us Express- EL388608922US 
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TECHNICAL FIELD 

The present invention relates to methods, processes and apparatus for measuring and anal- 
ysis of variables, provided that the definitions of the terms "variable" and "measuring" 
are adopted from the 7th edition of the International Patent Classification (IPC). This 
invention also relates to generic measurement systems and processes, that is, the proposed 
measuring arrangements are not specially adapted for any specific variables, or to one par- 
ticular environment. This invention also relates to methods and corresponding apparatus 
for measuring which extend to different applications and provide results other than instan- 
taneous values of variables. The invention further relates to post-processing analysis of 
measured variables and to statistical analysis. 
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BACKGROUND ART 

In a broad sense, the primary goal of a measurement can be defined as making a phe- 
nomenon available for human perception. Even when the results of measurements are used 
to automatically control machines and processes, the results of such control need to be 
5 meaningful, and thus the narrower technical meanings of a measurement still fall under 
the more general definition. In a technical sense, measurement often means finding a nu- 
j merical expression of the value of a variable in relation to a unit or datum or to another 

^ variable of the same nature. This is normally accomplished by practical implementation 
J of an idealized data acquisition system. The idealization is understood as a (simplified) 
j 10 model of such measuring process, which can be analyzed and comprehended by individuals. 

This analysis can either be performed directly through senses, or employ additional tools 
jj such as computers. When the measurement is reduced to a record, such record is normally 
3 expressed in discrete values in order to reduce the amount of information, and to enable 
3 storage of this record and its processing by digital machines. The reduction to a finite set 
"■ 15 of values is also essential for human comprehension. However, a physical embodiment of an 
idealized data acquisition system is usually an analog machine. That is, it is a machine with 
continuous action, where the components (mechanical apparatus, electrical circuits, optical 
devices, and so forth) respond to the input through the continuously changing parame- 
ters (displacements, angles of rotation, currents, voltages, and so forth). When the results 
20 of such implementation are reduced to numerical values, the uncertainties due to either 
limitations of the data acquisition techniques, or to the physical nature of the measured 
phenomenon, are often detached from these numerical values, or from each other. Ignoring 
the interdependence of different variables in the analyzed system, either intrinsic (due to 
their physical nature), or introduced by measuring equipment, can lead to misleading con- 
25 elusions. An example of such idealization of a measurement is its digital record, where the 
measurement is represented by a finite set of numbers. It needs to be pointed out that the 
digital nature of a record is preserved even if such record were made continuous in time, 
that is, available as a (finite) set of instantaneous values. 



2 



Atty Docket No. . 15252.2 
US Express: EL38860S922'uS 

Generally, measurement can be viewed as transformation of the input variable into 
another variable such that it can be eventually perceived, or utilized in some other man- 
ner. Measurement may consist of many intermediate steps, or stages between the incoming 
variable and the output of the acquisition system. For example, a TV broadcast can simplis- 
tically be viewed as (optical) measurement of the intensity of the light (incoming variable), 
where the output (image) is displayed on a TV screen. The same collectively would be true 
for a recorded TV program, although the intermediate steps of such a measurement will be 
different. 

Regardless of the physical nature of measuring processes, they all exhibit many common 
features. Namely, they all involve transformation and comparison of variables at any stage. 
Transformation may or may not involve conversion of the nature of signals (for instance, 
conversion of pressure variations into electric signals by a microphone in acoustic mea- 
surements), and transformation can be either hnear or nonhnear. Most transformations of 
variables in an acquisition system involve comparison as the basis for such transformations. 
Comparison can be made in relation to any external or internal reference, including the 
input variable itself. For example, simple linear filtering of a variable transforms the input 
variable into another variable, which is a weighted mean of the input variable either in time, 
space, or both. Here the comparison is made with the sample of the input variable, and 
the transformation satisfies a certain relation, that is, the output is the weighted average 
of this sample. An example of such filtering would be the computation of the Dow Jones 
Industrial Average. 

In measurements of discrete events, a particular nonlinear filtering technique stands out 
due to its important role in many applications. This technique uses the relative positions, 
or rank, of the data as a basis for transformation. For example, the salaries and the family 
incomes are commonly reported as percentiles such as current median salary for a certain 
profession. The rationale for reporting the median rather than the mean income can be 
illustrated as follows. Consider some residential neighborhood generating ten million dollars 
annually. Now, if someone from this neighborhood wins twenty millions in a lottery, this 
will triple the total as well as the mean income of the neighborhood. Thus reporting the 
mean family income will create an illusion of a significant increase in the wealth of individual 
families. The median income, however, .will remain unchanged and will reflect the economic 
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conditions of the neighborhood more accurately. As another simple example, consider the 
way in which a student's performance on a standardized test such as the SAT (Scholastic 
Aptitude Test) or GRE (Graduate Record Examination) is measured. The results are 
provided as a cumulative distribution function, that is, are quoted both as a "score" and 

5 as the percentile. The passing criterion would be the score for a certain percentile. This 
passing score can be viewed as the output of the "admission filter" . 

In digital signal processing, a similar filtering technique is commonly used and is re- 
ferred to as rank order or order statistic filtering. Unlike a smoothing filter which outputs 
a weighted mean of the elements in a sHding window, a rank order filter picks up an output 

10 according to the order statistic of elements in this window. See, for example, Arnold et al., 
1992, and Sarhan and Greenberg, 1962, for the definitions and theory of order statistics. 
Maximum, minimum, and median filters are some frequently used examples. Median filters 
are robust, and can remove impulse noise while preserving essential features. The discussion 
of this robustness and usefulness of median filters can be found in, for example, Arce et al., 

15 1986. These filters are widely used in many signal and image processing applications. See, 
for example, Bovik et al, 1983; Huang, 1981; Lee and Fam, 1987. Many examples can be 
found in fields such as seismic analysis Bednar, 1983, for example, biological signal process- 
ing Fiore et al., 1996, for example, medical imaging Ritenour et al., 1984, for example, or 
video processing Wischermann, 1991, for example. Maximum and minimum selections are 

20 also quite common in various applications Haralick et al., 1987, for example. 

Rank order filtering is only one of the applications of order statistic methods. In a 
simple definition, the phrase order statistic methods refers to methods for combining a 
large amount of data (such as the scores of the whole class on a homework) into a single 
number or small set of numbers that give an overall flavor of the data. See, for example, 

25 Nevzorov, 2001, for further discussion of different applications of order statistics. The main 
limitations of these methods arise from the explicitly discrete nature of their definition 
(see, for example, the definitions in Sarhan and Greenberg, 1962, and Nevzorov, 2001), 
which is in striking dissonance with the continuous nature of measurements. The discrete 
approach imposes the usage of algebraic rather than geometric tools in order statistics, and 

30 thus limits both the perception of the results through the geometric interpretation and the 
applicability of differential methods of analysis. 
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Order statistics of a sample of a variable is most naturally defined in terms of the 
cumulative distribution function of the elements composing this sample see David, 1970, 
for example, which is a monotonic function. Thus computation of an order statistic should 
be equivalent to a simple task of finding a root of a monotonic function. However, the 
5 cumulative distribution of a discrete set is a discontinuous function, since it is composed of 
a finite number of step functions (see Scott, 1992, for example). As a result, its derivative 
(the density function) is singular, that is, composed of a finite number of impulse functions 
such as Dirac (5-function (see, for example, Dirac, 1958, p. 58-61, or Davydov, 1988, p. 609- 
612, for the definition and properties of the Dirac 5-function). When implementing rank 
10 order methods in software, this discontinuity of the distribution function prevents us from 
using efficient methods of root finding involving derivatives, such as the iNewton-Raphson 
method (see Press et al., 1992, and the references therein for a discussion of root finding 
methods). In hardware, the inabihty to evaluate the derivatives of the distribution function 
disallows analog implementation. Even though for a contmuous-Ume signal the distribution 
15 function may be continuous in special cases (since now it is an average of an infinitely large 
number of step functions), the density function is still only piecewise continuous, since 
every extremum in the sample produces singularity in the density function (Nikitin, 1998, 
Chapter 4, for example). In fact, the nature of a continuous-time signal is still discrete, since 
its instantaneous and even time averaged densities are still represented by impulse functions 
20 (Nikitin, 1998, for example). Thus the time continuity of a signal does not automatically 
lead to the continuity of the distribution and the density functions of a sample of this signal. 

Following from their discrete nature, the limitations of the existing rank order methods 
(rank order filtering as well as other methods based on order statistics) can roughly be 
divided into two categories. The first category deals with the issues of the implementation 
25 of these methods, and the second one addresses the limitations in the applicability. The 
implementation of the order statistics methods can in turn be divided into two groups. 
The first group realizes these methods in software on sequential or parallel computers (see 
Juhola et al., 1991, for example). The second one implements them on hardware such as 
Very Large Scale Integration (VLSI) circuits (see Murthy and Swamy, 1992, for example). 
30 In software implementation, the basic procedure for order statistics calculation is com- 
parison and sorting. Since sorting can be constructed by selection, which is an operation 
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linear in complexity, the algorithms for finding only a specific rank (such as median) are 
more effective than the algorithms for computation of arbitrary statistics (Pasian, 1988, 
for example). In addition, the performance of rank order calculations can be improved by 
taking advantage of the running window where only a minor portion of the elements are 
deleted and replaced by the same number of new elements (Astola and Campbell, 1989, 
for example). Regardless of the efficiency of particular algorithms, however, all of them 
quickly become impractical when the size of the sample grows, due to the increase in both 
computational intensity and memory requirements. 

The hardware implementation of rank order processing has several main approaches, 
such as systolic algorithms (Fisher, 1984, for example), sorting networks (Shi and Ward, 
1993, and Opris, 1996, for example), and radix (binary partition) methods (Lee and Jen, 
1993, for example). The various hardware embodiments of the order statistics methods, 
however, do not overcome the intrinsic limitations of the digital approach arising from the 
discontinuous nature of the distribution function, such as inefficient rank finding, difficulties 
with processing large samples of data, and inability to fully explore differential techniques 
of analysis. It needs to be pointed out that the differential methods allow studying the 
properties "at a point" , that is, the properties which depend on an arbitrary small neigh- 
borhood of the point rather than on a total set of the discrete data. This offers more 
effective technical solutions. Several so-called "analog" solutions to order statistic filter- 
ing have been proposed see Jarske and Vainio, 1993, for example, where the term "analog" 
refers to the continuous (as opposed to quantized) amplitude values, while the time remains 
discrete. Although a definition of the continuous-time analog median filter has been known 
since the 1980's (see Fitch et al., 1986), no electronic implementations of this filter have 
been introduced. Perhaps the closest approximation of the continuous-time analog median 
filter known to us is the Unear median hybrid (LMH) filter with active RC hnear subfilters 
and a diode network (Jarske and Vainio, 1993, for example). 

The singular nature of the density functions of discrete variables does not only impede 
both software and hardware implementations of rank order methods, but also constrains the 
applicability of these methods (for example, their geometric extension) to signal analysis. 
The origin of these constraints lies in the contrast between the discrete and the contin- 
uous: "The mathematical model of a separate object is the unit, and the mathematical 
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model of a collection of discrete objects is a sura of units, which is, so to speak, the image 
of pure discreteness, purified of all other qualities. On the other hand, the fundamental, 
original mathematical model of continuity is the geometric figure; ..." (Aleksandrov et al., 
1999, V. I, p. 32). Even simple time continuity of the incoming variable enables differenti- 

5 ation with respect to time, and thus expands such applicability to studying distributions 
of local extrema and crossing rates of signals (Nikitin et al., 1998, for example), which 
can be extremely useful characteristics of a dynamic system. However, these distributions 
are still discontinuous (singular) with respect to the displacement coordinates (thresholds). 
Normally, this discontinuity does not restrain us from computing certain integral charac- 

10 teristics of these distributions, such as their different moments. However, many useful tools 
otherwise applicable to characterization of distributions and densities are unavailable. For 
instance, in studies of experimentally acquired distributions the standard and absolute devi- 
ations are not reliable indicators of the overall widths of density functions, especially when 
these densities are multimodal, or the data contain so-called outliers. A well-known quan- 

15 tity Full Width at Half Maximum (FWHM) (e.g., Zaidel' et al, 1976, p.l8), can characterize 
the width of a distribution much more reliably, even when neither standard nor absolute 
deviation exists. The definition of FWHM, however, requires that the density function be 
continuous and finite. One can introduce a variety of other useful characteristics of distri- 
butions and density functions with clearly identifiable geometrical and physical meaning, 

20 which would be unavailable for a singular density function. An additional example would 
be an a-level contour surface (Scott, 1992, p. 22), which requires both the continuity and 
the existence of the maximum or modal value of the density function. 

Discontinuity of the data (and thus singularity of density functions) is not a direct re- 
sult of measurements but rather an artifact of idealization of the measurements, and thus a 

25 digital record should be treated simply as a sample of a continuous variable. For example, 
the threshold discontinuity of digital data can be handled by convolution of the density 
function of the discrete sample with a continuous kernel. Such approximation of the "true" 
density is well known as Kernel Density Estimates (KDE) (Silverman, 1986, for example), 
or the Parzen method (Parzen, 1967, for example). This method effectively transforms a 

30 digital set into a threshold continuous function and allows successful inference of "true" 
distributions from observed samples. See Lucy, 1974, for the example of the rectification 



7 



Atty. Docket No. . 16252.2 
First Named Inventor. Nikitin 
US Express EL388608922US 

of observed distributions in statistical astronomy. The main limitation of the KDE is that 
the method primarily deals with samples of finite size and does not allow treatment of 
spatially and temporally continuous data. For example, KDE does not address the time 
dependent issues such as order statistic filtering, and does not allow extension of the contin- 

5 uous density analysis to intrinsically time dependent quantities such as counting densities. 
Another important limitation of KDE is that it fails to recognize the importance of and 
to utilize the cumulative distribution function for analysis of multidimensional variables. 
According to David W. Scott (Scott, 1992, page 35), "...The multivariate distribution 
function is of little interest for either graphical or data analytical purposes. Furthermore, 

10 ubiquitous multivariate statistical apphcations such as regression and classification rely on 
direct manipulation of the density function and not the distribution function" . Some other 
weaknesses of KDE with respect to the present invention will become apparent from the 
further disclosure. 

Threshold, spatial, and temporal continuity are closely related to our inability to conduct 

15 exact measurements, for a variety of reasons ranging from random noise and fluctuations to 
the Heisenberg uncertainty. Sometimes the exact measurements are unavailable even when 
the measured quantities are discrete. An example can be the "pregnant chad" problem in 
counting election votes. As another example, consider the measurement of the energy of a 
charged particle. Such measurement is normally carried out by means of discriminators. 

20 An ideal discriminator will register only particles with energies larger than its threshold. 
In reality, however, a discriminator will register particles with smaller energies as well, and 
will not detect some of the particles with larger energies. Thus there will be uncertainties 
in our measurements. Such uncertainties can be expressed in terms of the response function 
of the discriminator. Then the results of our measurements can be expressed through the 

25 convolution of the "ideal" measurements with the response function of the discriminator 
(Nikitin, 1998, Chapter 7, for example). Even for a monoenergetic particle beam, our 
measurements will be represented by a continuous curve. Since deconvolution is at least an 
impractical, if not impossible, way of restoring the "original" signal, the numerical value 
for the energy of the incoming particles will be deduced from the measured density curve 

30 as, for example, its first moment (ZaideP et al., 1976, pp. 11-24, for example). 

A methodological basis for treatment of an incoming variable in terms of its continuous 
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densities can be found in fields where the measurements are taken by an analog action 
machine, that is, by a probe with continuous (spatial as well as temporal) impulse response, 
such as optical spectroscopy (see Zaidel' et al, 1976, for example). The output of such a 
measuring system is described by the convolution of the impulse response of the probe with 
5 the incoming signal, and is continuous even for a discrete incoming signal. For instance, 
the position of a certain spectral line measured by a monochromator is represented by a 
smooth curve rather than by a number. If the reduction of the line's position to a number is 
needed, this reduction is usually done by replacing the density curve by its modal, median, 
or average value. 

10 The measurement of variables and analysis of signals often go hand-in-hand, and the 

distinction between the two is sometimes minimal and normally well understood from the 
context. One needs to understand, however, that a "signal", commonly, is already a result 
of a measurement. That is, a "signal" is already a result of a transformation (by an acqui- 
sition system) of one or many variables into another variable (electrical, optical, acoustic, 

15 chemical, tactile, etc.) for some purpose, such as further analysis, transmission, directing, 
warning, indicating, etc. The relationship between a variable and a signal can be of a 
simple type, such that an instantaneous value of the variable can be readily deduced from 
the signal. Commonly, however, this relationship is less easily decipherable. For example, 
a signal from a charged particle detector is influenced by both the energies of the particles 

20 and the times of their arrival at the sensor. In order to discriminate between these two 
variables, one either needs to use an additional detector (or change the acquisition param- 
eters of the detector), or to employ additional transformation (such as differentiation) of 
the acquired signal. The analysis of the signal is thus a means for gathering information 
about the variables generating this signal, and ultimately making the phenomenon available 

25 for perception, which is the goal. The artificial division of this integral process into the 
acquisition and the analysis parts can be a serious obstacle in achieving this goal. 

In the existing art, the measurement is understood as reduction to numbers, and such 
reduction normally takes place before the analysis. Such premature digitization often un- 
necessary complicates the analysis. The very essence of the above discussion can be revealed 

30 by the old joke that it might be hard to divide three potatoes between two children un- 
less you make mashed potatoes. Thus we recognize that the nature of the difficulties with 
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implementation and applicability of order statistics methods in analysis of variables lies in 
the digital approach to the problem. By digitizing, we lose continuity. Continuity does not 
only naturally occur in measurements conducted by analog machines, or arise from consid- 
eration of uncertainty of measurements. It is also important for perception and analysis 

5 of the results of complex measurements, and essential for geometrical and physical inter- 
pretation of the observed phenomena. The geometric representation makes many facts of 
analysis "intuitive" by analogy with the ordinary space. By losing continuity, we also lose 
differentiability, which is an indispensable analytical tool since it allows us to set up differ- 
ential equations describing the studied system: ". . . In order to determine the function that 

10 represents a given physical process, we try first of all to set up an equation that connects 
this function in some definite way with its derivatives of various orders" (Aleksandrov et al., 
1999, V. I, p. 119). 

The origin of the limitations of the existing art can thus be identified as relying on 
the digital record in the analysis of the measurements, which impedes the geometrical 
15 interpretation of the measurements and leads to usage of algebraic rather than differential 
means of analysis. 
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DISCLOSURE OF INVENTION 



BRIEF SUMMARY OF THE INVENTION 



As was stated in the description of the existing art, the digital approach hmits both the 
5 geometrical interpretation of the measurements and prevents usage of the differential means 
of analysis. In this invention, we present an analog solution to what is usually handled by 
digital means. We overcome the deficiencies of the prior approach by considering, instead 
of the values of the variables, such geometrical objects as the densities of these variables in 
their threshold space. The applicability of the differential analysis is achieved by either (1) 
10 preserving, whenever possible, the original continuity of the measurement in the analysis, 
or (2) restoring continuity of discrete data through convolution with a continuous kernel 
which represents the essential qualities of the measuring apparatus. Our approach offers 
(a) improved perception of the measurements through geometrical analogies, (b) effective 
solutions of the existing computational problems of the order statistic methods, and (c) 
15 extended applicability of these methods to analysis of variables. In the subsequent disclosure 
we will demonstrate the particular advantages of the invention with respect to the known 
art. 

In this invention, we address the problem of measuring and analysis of variables, on a 
consistent general basis, by introducing threshold densities of these variables. As will be 

20 described further in detail, the threshold densities result from averaging of instantaneous 
densities with respect to thresholds, space, and time. Since this averaging is performed by a 
continuous kernel (test function), it can be interpreted as analog representation, and thus 
we adopt the collective designation Analysis of Variables Through Analog Representation 
(AVATAR) for various components of the invention. The interdependence of the variables 

25 in the measurement system is addressed by introducing the modulated threshold densities 
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of these variables, which result from the consideration of the joint densities of the inter- 
dependent variables in their combined threshold space. The particular way in which these 
densities are introduced leads to the densities being continuous in thresholds, space, and 
time, even if the incoming variables are of discrete nature. This approach allows us to 

5 successfully address the limitations of the prior art identified earlier, opens up many new 
opportunities for expanding the applicability of rank order analysis of variables, and pro- 
vides a means for efficient implementation of this analysis in both hardware and software. 

In order to convey the inventive ideas clearly, we adopt the simplified model for mea- 
surements as follows (see the analogous ideal system in Nikitin et al., 1998, for example). 

10 A variable is described in terms of displacement coordinates, or thresholds, as well as in 
terms of some other coordinates such as spatial coordinates and physical time. The values 
of these coordinates are measured by means of discriminators and/or differential discrim- 
inators (probes). An ideal discriminator with threshold D returns the value "1" if the 
measured coordinate exceeds D, "1/2" if it is equal to D, and it returns zero otherwise. 

15 Thus the mathematical expression for an ideal discriminator is the Heaviside unit step func- 
tion of the difference between the threshold and the input coordinate (see Nikitin et al., 
1998, and Nikitin, 1998, for example). Although ideal discriminators can be a useful tool 
for analysis of a measurement process, different functions of • thresholds need to be em- 
ployed to reflect the vagaries of real measurements (see discussion in Nikitin et al., 1998, 

20 for example). In this invention, the peculiarities of "real" discriminators are reflected by 
introducing uncertainty into the mathematical description of the discriminator, such that 
the returned value is in the range zero to one, depending on the values of the input co- 
ordinate. As described further in this disclosure, the introduction of such uncertainty can 
be interpreted as averaging with respect to a threshold test function, or threshold averag- 

25 ing. When the input-output characteristic of the discriminator is a continuous function, 
then differentiability of the output with respect to threshold is enabled. In addition, if the 
characteristic of the discriminator is a monotonic function, the rank relations of the input 
signal are preserved. If the original input signal is not differentiable with respect to time 
(e.g., the input signal is discontinuous in time), differentiability with respect to time can 

30 always be enabled by introducing time averaging into the acquisition system, where under 
time averaging we understand a suitable convolution transform with a continuous-time ker- 
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nel. Likewise, differentiability with respect to spatial coordinates can be enabled by spatial 
averaging. 

The mathematical expression for the response of an ideal probe, or differential discrim- 
inator, is the Dirac 5-function of the difference between the displacement and the input 
variable, which is the instantaneous density of the input variable. As follows from the 
properties of the Dirac (5-function (see, for example, Dirac, 1958, p. 58-61, and Davydov, 
1988, p. 609-612, for the definition and properties of the Dirac 5-function), the output of 
the "real" probe is thus the convolution of the instantaneous density with the input- output 
characteristic of the differential discriminator, which is equivalent to the threshold averag- 
ing of the instantaneous density. When this output is subsequently averaged with respect 
to space and time, the result is the Threshold- Space- Time Averaged Density. 

Notice that the transition from the ideal to real probes and discriminators preserves 
the interpretation of their responses as the threshold density and cumulative distribution, 
respectively. For example, the spectrum acquired by an optical spectrograph can be con- 
sidered the energy density regardless the width of the slits of its monochromator. Thus a 
particular way of utilizing the discriminators and probes in this invention is essentially a 
method of transformation of discrete or continuous variables, and/or ensembles of variables 
into normalized continuous scalar fields, that is, into objects with mathematical properties 
of density and cumulative distribution functions. In addition to dependence on the displace- 
ment coordinates (thresholds), however, these objects can also depend on other parameters, 
including spatial coordinates (e.g., ii the incoming variables are themselves scalar or vector 
fields), and/or time (if the variables depend on time). For the purpose of this disclosure, 
the terms "space" and "time" are used to cover considerably more than their primary or 
basic meaning. "Time" should be understood as a monotonic scalar, continuous or dis- 
crete, common to all other analyzed variables, which can be used for sequential ordering 
of the measurements. "Space" is thus all the remaining coordinates which are employed 
(as opposed to sufficient) to govern the values of the input variables. It is important to 
note that the use of the invented transformation makes all existing techniques for statistical 
processing of data, and for studying probability fluxes, applicable to analysis of these vari- 
ables. Moreover, the transformations of the measured variables can be implemented with 
respect to any reference variable, or ensemble of reference variables. In this disclosure, we 
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consider two basic transformations with respect to the reference variable, which we refer 
to as normalization and modulation. The definitions of these transformations will be given 
later in the disclosure. In both of these transformations, the behavior of the input variable 
is represented in terms of behavior (and units) of the reference variable. Thus, the values 
of the reference variable provide a common unit, or standard, for measurement and com- 
parison of variables of different nature, for assessment of the mutual dependence of these 
variables, and for evaluation of the changes in the variables and their dependence with time. 
For example, dependence of economic indicators on social indicators, and vice versa, can 
be analyzed, and the historical changes in this dependence can be monitored. When the 
reference variable is related in a definite way to the input variable itself, these additional 
transformations (that is, normalization and modulation) provide a tool for analysis of the 
interdependence of various properties of the input variable. 

Among various embodiments of the invention, several are of particular importance for 
analysis of variables. These are the abihty to measure (or compute from digital data) (1) 
quantile density, (2) quantile domain, and (3) quantile volume for a variable. Quantile 
density indicates the value of the density likely to be exceeded, quantile domain contains 
the regions of the highest density, and quantile volume gives the (total) volume of the 
quantile domain. The definitions of these quantities and a means of their implementation 
are unavailable in the existing art. Detailed definitions and description of these quantities 
will be given later in the disclosure. As another consequence of the proposed transformation 
of variables into density functions, the invention enables measurements of currents, or fluxes 
of these densities, providing a valuable tool for analysis of the variables. 

Another important embodiment of AVATAR is rank normalization of a variable with 
respect to a cumulative distribution function, generated by another variable or ensemble of 
variables. The rank normalization of variables can be used for processing and analysis of 
different time ordered series, ensembles, scalar or vector fields, and time independent sets 
of variables, especially when the investigated characteristics of the variable are invariant 
to a monotonic transformation of its values, that is, to a monotonic transformation of the 
thresholds. This normalization can be performed with respect to a reference distribution of 
an arbitrary origin, such as the distribution provided by an external reference variable, or 
by the input variable itself For example, the reference variable can be a random process 
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with the parameters determined by the input variable. In this case, the reference variable 
provides a "container" in the threshold space where the input variable is likely to be found. 
More importantly, the rank normalization allows computation or measurement of integral 
estimators of differences between two distributions (or densities) as simple time and/or 
space averages. 

Another important usage of the rank normalization is as part of preprocessing of the 
input variable, where under preprocessing we understand a series of steps (e.g., smoothing) 
in the analysis prior to applying other transformations. Since in AVATAR the extent of the 
threshold space is determined by the reference variable, the rank normalization allows as 
to adjust the resolution of the acquisition system according to the changes in the threshold 
space, as the reference variable changes in time. Such adjustment of the resolution is the 
key to a high precision of analog processing. 

While rank normalization reflects the rank relations between the input and the reference 
variables, the modulated threshold density describes the input variable in terms of the rate 
of change of the reference variable at a certain threshold. As will be clarified further 
in the disclosure, the modulated threshold densities arise from the consideration of the 
joint densities of the input and the reference variables in their combined threshold space. 
Instead of analyzing such joint variable in its total threshold space, however, we consider 
the behavior of the input variable in the threshold space of the reference variable only. 
The modulated densities allow us to investigate the interdependence of the input and the 
reference variable by comparison of the time averages of the input variable at the thresholds 
of the reference variable with the simple time average of the input variable. As has been 
mentioned earlier, the important special case of modulated densities arises when there is a 
definite relation between the input and the reference variables. In this disclosure, we will 
primarily focus on the densities arising from the two particular instances of this relation, 
to which we will further refer as the amplitude and the counting densities. Since the 
definiteness in the relation eliminates the distinction between the input and the reference 
variables, such special cases of the modulated densities will be regarded simply as modulated 
threshold densities of the input variable. 

The invention also allows the transformation which can be interpreted as rank filtering 
of variables. That is, it enables the transformation of a variable into another variable, the 
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value of which at any given space and time is a certain quantile of a modulated cumulative 
distribution function generated by the input and reference variables. Thus in AVATAR 
the output of such a filter has simple and clear interpretation as a level line of such a 
distribution in the time-threshold plane. One needs to notice that the output of a rank 
filter as defined in the existing digital signal processing methods, will correspond to the 
discrete points on a level line drawn for an amplitude distribution only. Thus the rank 
filtering defined in AVATAR extends beyond the known applicability of rank filtering. It is 
also important that, in this invention, such filtering process is implemented by differential 
means and thus conducted without sorting. The invention also provides a means for finding 
(selecting) the rank of a time dependent or static variable or ensemble of variables without 
sorting. Such analog rank selection permits analog emulation of digital rank filters in an 
arbitrary window, of either finite or infinite type, with any degree of precision. Moreover, 
the rank selection is defined for modulated distributions as well, which extends the applica- 
bility of rank filtering. The continuous nature of the definitions and the absence of sorting 
allows easy implementation and incorporation of rank selecting and filtering in analog de- 
vices. Rank filtering and selecting also provide alternative embodiments for comparison of 
variables with respect to a common reference variable, and for detection and quantification 
of changes in variables. 

Based on the embodiments of AVATAR discussed above, we can define and implement 
a variety of new techniques for comparison of variables and for quantification of changes 
in variables. By determining distributions from the signal itself, and/or by providing a 
common reference system for variables of different natures, the invention provides a robust 
and efficiently applied solution to the problem of comparing variables. It is important to 
note that the invention enables comparison of one-dimensional as well as multivariate den- 
sities and distributions by simple analog machines rather than through extensive numerical 
computations. 

The particular advantages of AVATAR stem from the fact that the invention is based 
on the consideration of real acquisition systems while most of previous techniques assume 
idealized measurement processes. Even if the incoming variable is already a digital record 
(such as the result of an idealized measurement process), we restore this record to a con- 
tinuous form by a convolution with a probe continuous in thresholds, space, and time. It is 
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important to realize that the invention is also applicable to measurements of discrete data. 
The main distinction between such measurements and the restoration of a digital record 
lies in the realization that the former are normally taken by continuous action machines, 
and thus the results of such measurements need to be reformulated, that is, "mapped" into 
analog domain. As an example, consider the task of measuring the amphtude (energy) 
distribution of a train of charged particles. Since this measurement is usually taken by 
means of sensors with finite time response, the problem of measuring the amplitude den- 
sity is naturally restated as the problem of finding the distribution of local extrema of a 
continuous-time signal (Nikitin et al., 1998). This distribution can in turn be found through 
the derivative of the crossing rates with respect to threshold (Nikitin et al., 1998). 

The invention can be implemented in hardware devices as well as in computer codes 
(software). The applications of the invention include, but are not limited to, analysis of a 
large variety of technical, social, and biologic measurements, traffic analysis and control, 
speech and pattern recognition, image processing and analysis, agriculture, and telecom- 
munications. Both digital and analog implementations of the methods can be used in 
various systems for data acquisition and analysis. All the above techniques, processes and 
apparatus are applicable to analysis of continuous (analog) variables as well as discrete 
(digital). By analyzing the variables through their continuous distributions and the density 
functions, the invention overcomes the limitations of the prior art by (a) improving percep- 
tion of the measurements through geometrical analogies, (b) providing effective solutions 
to the existing computational problems of the order statistic methods, and (c) extending 
the applicability of these methods to analysis of variables. 

Further scope of the applicability of the invention will be clarified through the detailed 
description given hereinafter. It should be understood, however, that the specific examples, 
while indicating preferred embodiments of the invention, are presented for illustration only 
Various changes and modifications within the spirit and scope of the invention should 
become apparent to those skilled in the art from this detailed description. Furthermore, 
all the mathematical expressions and the examples of hardware implementations are used 
only as a descriptive language to convey the inventive ideas clearly, and are not limitative 
of the claimed invention. 
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Terms and Definitions 
WITH Illustrative Examples 

For convenience, the essential terms used in the subsequent detailed description of the 
invention are provided below. These terms are listed along with their definitions adopted 
for the purpose of this disclosure. Examples clarifying and illustrating the meaning of 
5 the definitions are also provided. Note that the equations in this section are numbered 
separately from the rest of this disclosure. 

1 Variable 

For the purpose of this disclosure, we define a variable as an entity x which, when ex- 
pressed in figures (numbers), can be represented by one of the mathematical expressions 

10 listed below. Throughout this disclosure, the standard mathematical terms such as vector, 
scalar, or field mostly preserve their commonly acceptable mathematical and/or physical 
interpretation. The specific meaning of most of the common terms will be clarified through 
their usage in the subsequent detailed description of the invention. Notice that the most 
general form of a variable adopted in this disclosure is an Ensemble of Vector Fields. All 

15 other types of variables in this disclosure can be expressed through various simplifications of 
this general representation. For example, setting the ensemble weight n(jj,) in the expression 
for an ensemble of vector fields to the Dirac 5-function d(^) reduces said expression to a 
single Vector Field variable. Further, by eliminating the dependence of the latter on spatial 
coordinates (that is, by setting the position vector a = constant), said single vector field 

20 variable reduces to a single Vector variable. Notice also that while ensembles of variables 
are expressed as integrals/sums of the components of an ensemble, it should be understood 
that individual components of an ensemble are separately available for analysis. 

1. Single Variable. A single variable in this disclosure can be a vector/scalar variable, 
or a vector/scalar field variable. 
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(a) A Vector Field variable can be expressed as 

x = x(a,t), (D-i) 
where a is the vector of spatial coordinates, and t is the time coordinate. 

Representative Examples: (1) A truecolor image can be expressed by a vector field 
5 X = x(a,t), where the color is described by its coordinates in the three-dimensional color 
space (red, green, and blue) at the position a. (2) The combination of both the color 
intensity of a monochrome image and the rate of change of said intensity can be expressed 
by a vector field x = (a; {a, t), i;(a, t)). 

(b) A Scalar Field variable can be expressed as 

10 

x = a;(a,i), (D-2) 
where a is the vector of spatial coordinates, and t is the time coordinate. 

Representative Example: A monochrome image at a given time is determined by the 
intensity of the color at location a, and thus it is conveniently described by a scalar field 
X = a;(a, i). 

15 (c) A Vector variable can be represented by the expression 

X = x(i) , (D-3) 
where t is the time coordinate. Notice that the components of a vector variable do not have 
to be of the same units or of the same physical nature. 

Representative Examples: (1) The position of a vehicle in a traflic control problem can 
20 be described by a single vector variable x = x(t). (2) The position and the velocity of a 
vehicle together can be described by a single vector variable x = (x(i),x(^)). 
(d) A Scalar variable can be represented by the expression 

X = , (D-4) 

where t is the time coordinate. 
25 Representative Example: The current through an element of an electrical circuit can be 

expressed as a single scalar variable x = 

2. Ensemble of Variables, Several different variables of the same nature can be con- 
sidered as a (composite) single entity designated as an ensemble of variables. An individual 
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variable in an ensemble is a component of the ensemble. The relative contribution of a 
component in the ensemble is quantified by a weight n{ii) of the component. 

(a) Ensemble of Vector Fields: 

^ " I^J"^ ^^^^ ^) ' (D-5) 
where n{fj) d^i is the wez^/^i of the /ith component of the ensemble such that /f^d/in(/i) = 
A^, a is the vector of spatial coordinates, and t is the time coordinate. 

Representative Example: A truecolor image can be expressed by a vector field x = 
x(a, t), where the color is described by its coordinates in the three-dimensional color space 
(red, green, and blue), at the position a. A "compound" image consisting of a finite or 
infinite set of such truecolor images, weighted by the weights can be viewed as an 

ensemble of vector fields. For example, such a compound image can be thought of as a 
statistical average of the video recordings taken by several different cameras. 

(b) Ensemble of Scalar Fields: 

^ = f_jLIJ^n{j2)x^{a,t) , (D-6) 
where n{^) dfi is the weight of the /^th component of the ensemble such that /^d//n(/i) = 
iV, a is the vector spatial coordinates, and t is the time coordinate. 

Representative Example: A monochrome image at a given time is determined by the 
intensity of the color at location a, and thus it is conveniently described by a scalar field 
X = x{a.,t). A "compound" image consisting of a finite or infinite set of such monochrome 
images, weighted by the weights n(/:x), can be viewed as an ensemble of scalar fields. For 
example, such a compound image can be thought of as a statistical average of the video 
recordings taken by several different cameras. 

(c) Ensemble of Vectors: 

X = /^d/i K^{t) , (D-7) 

where n(/i) d^ is the weight of the juth component of the ensemble such that /^d^n(/x) = 
A^, and t is the time coordinate. 

Representative Example: A variable expressing the average position of different ve- 
hicles in a traffic control problem can be described by an ensemble of vector variables 
X = /^d/in(//) x^(i), where n{fi) = n(/i) E^Ii S{fx-i). 
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(d) Ensemble of Scalars: 

x = J°°diun{fj.)x^{t) , (D-8) 

where n(/i) d/j, is the weight of the {ith component of the ensemble such that /f^dyU.n(yu) = 
iV, and t is the time coordinate. 
5 Representative Example: The total current through N elements of an electrical circuit 

can be expressed as an ensemble of scalar variables x = /f^d/i7z(/i) where n{iJ,) = 

Notice that the most general form of a variable among those listed above is the Ensemble 
of Vector Fields. All other types of variables in this disclosure can be expressed through var- 

10 ious simplifications of this general representation. For example, setting the ensemble weight 
n(fj,) in the expression for an ensemble of vector fields to the Dirac (5-function reduces 
said expression to a single Vector Field variable. Further, by eliminating the dependence 
of the latter on spatial coordinates (that is, by setting the position vector a = constant), 
said single vector field variable reduces to a single Vector variable. Notice also that while 

15 ensembles of variables are expressed as integrals/sums of the components of an ensemble, it 
should be understood that individual components of an ensemble are separately available 
for analysis. 



2 Threshold Filter 

In this disclosure, we define a Threshold Filter as a continuous action machine (a physical 
20 device, mathematical function, or a computer program) which can operate on the differ- 
ence between a Displacement Variable D and the input variable x, and the result of such 
operation can be expressed as a scalar function of THE Displacement Variable, that is, as 
a value at a given displacement D. The dependence of the output of a Threshold Filter on 
the input is equivalent to those of a probe (smoothing threshold filter) or a discriminator 
25 (integrating threshold filter) as specified below. 

1. Probe. A continuous action machine (a physical device, mathematical function, or 
a computer program) which can operate on the difference between a Displacement Variable 
D and the input variable x, and the result of such operation can be expressed as a scalar 
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function of the Displacement Variable, that is, as a value at a given displacement D. The 
dependence of the output y of a probe on its input r is equivalent to the following expression: 

0 < y = /R(r) , (D-9) 
where y is a scalar, r is a vector or a scalar, R is a scalar width parameter of the probe, 
and the test function /r is a continuous function satisfying the conditions 
i /'Vr/R(r) = l, and 

< ^-"^ (D-10) 

where 5{r) is the Dirac (5-function. 

Representative Examples: (1) The Gaussian test function 

Mr) = nS..^..,W^ exp[-|:(^)^] (D-U) 

can act as a probe. In Eq. (D-11), the response of the probe /R(r) to the vector input 

r= (ri,...,r„) is a product of the responses of the probes So^J^A A (^i) = ad to 
the components rj of the input vector, (2) Fig. 6 illustrates an optical threshold smoothing 
filter (probe). This probe consists of a point light source S and a thin lens with the focal 
length /. The lens is combined with a gray optical filter with transparency described by 
/2/(x). Both the lens and the filter are placed in an XOY plane at a distance 2/ from the 
source S. The lens-filter combination can be moved in the XOY plane by the incoming 
signal r so that the center of the combination is located at in this plane. Then the 
output of the filter is proportional to the intensity of the light measured at the location 
D = (LI J,, Dy) in the D^-O-Dy plane parallel to the XOY plane and located at the distance 
R from the image S' of the source S (toward the source). That is, the output of this filter 
can be described by /r(D - r). 

2. Discriminator. A continuous action machine (a physical device, mathematical 
function, or a computer program) which can operate on the difference between a Dis- 
placement Variable D and the input variable x, and the result of such operation can be 
expressed as a scalar function of the Displacement Variable, that is, as a value at a given 
displacement D. The dependence of the output y of a discriminator on the input r is 
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obtained by integrating the response of the respective probe, that is, it is related to the 
input-output characteristic of the respective probe as 

0 < y = ^r(x) = £j^r Mr) < 1 . (D-12) 

5 Representative Examples: (1) The integral of a Gaussian test function 

•^R W = £^d"r Ur) = n -^AA (x,) = 2- fl erfc (^) , (D-13) 

where erfc(a;) is the complementary error function, can act as a discriminator. In Eq. (D- 
13), the response of the discriminator J^k{^) to the vector input x = {xi,...,Xn) is a 
product of the responses of the discriminators J^AoAxi) = |erfc (^) to the components 
10 Xi of the input vector. (2) By replacing the transparency function /2/(x) of the gray filter 
with J^2f{x), the optical probe shown in Fig. 6 is converted into a discriminator with the 
output J^jiCD - r). 

3 Displacement Variable 

A Displacement Variable is the argument of a function describing the output of a Threshold 
15 Filter. For example, if the Threshold Filter is an amplifier operating on the difference 
between two electrical signals D and x{t), and an output of the amphfier is described as a 
function of the input signal D, this signal D is a Displacement Variable. 

4 Modulating Variable 

A Modulating Variable is a unipolar scalar field variable K = K{a, t) which can be appHed 
20 to the output of a Threshold Filter in a manner equivalent to multiplication of said output 
by the Modulating Variable. 
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Representative Example: Imagine that the point light source S in Fig. 6 is produced 
by an incandescent lamp powered by a unit current. If we now power the lamp by the 
current K{t), the output of the threshold filter will be modulated by the Modulating Vari- 
able \K{t)\. 

5 Averaging Filter 

5 

An Averaging Fitter is a continuous action machine (a physical device, mathematical 
function, or a computer program) which can operate on a variable x(a,^), and the result 
of such operation can be expressed as convolution with a test function /R(a) and a time 
impulse response function h(t; T), namely as 

10 (x(r,s))Ji = ((x(r,.)>^)^ = ((x(r,.))^)^ = 

where a is the position vector (vector of spatial coordinates), and t is the time coordinate. 
Thus an averaging filter performs both spatial and time averaging. We shall call the product 
H^; T) /R(a) the impulse response of the Averaging Filter. 

1. Time Averaging Filter. An averaging filter which performs only time averaging 
15 is obtained by setting the spatial impulse response (test function) of the averaging filter 
to be equal to the Dirac 5-function, /R(a) = 5(a). The result of the operation of a time 
averaging filter with an impulse response h{t;T) on a variable x(a, i) can be expressed by 
the convolution integral 

(x(a, 3)>J = ^"^ds h{t - s; T) x(a, s) , (D-15) 
20 where T is the width (time scale) parameter of the filter. For two filters with the width 
parameters T and i\T such that T » AT, the former filter is designated as wide, and the 
latter as narrow. 

Representative Example: An image formed on a luminescent screen coated with lu- 
minophor with the afterglow half-time Ty^ = T In (2) is time averaged by an exponentially 
25 forgetting Time Averaging Filter h{t] T) = e~^ e{t) /T. 
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2. Spatial Averaging Filter. An averaging filter which performs only spatial averaging 
is obtained by setting the time impulse response of the averaging filter to be equal to the 
Dirac 5-function, h{t-T) = 5{t). The result of the operation of a spatial averaging filter 
with an impulse response (test function) /R(a) on a variable x(a,t) can be expressed by 
the convolution integral 

(x(r, t)Y^ = j^JJ'x f^{aL - r) x(r, t) , (D-16) 
where R is the width parameter of the filter. For two filters with the width parameters 
R and AR such that R > AR, the former filter is designated as wide, and the latter as 
narrow . 

Representative Example: A monochrome image given by the matrix Z = Zij(t) 
can be spatially averaged by a smoothing filter Wmn, T.mn'^mn = 1 as 

Some other terms and their definitions which appear in this disclosure will be provided 
in the detailed description of the invention. 
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NOTATIONS 

For convenience, lists of the acronyms and selected notations used in tlie detailed 
description of the invention are provided below. 



Selected Acronyms and Where They First Appear 



5 

AVATAR Analysis of Variables Through Analog Representation, page 11 

MTD Modulated Threshold Density, page 51 

MRT Mean at Reference Threshold, page 52 

MCTD Modulated Cumulative Threshold Distribution, page 66 

10 ARN Analog Rank Normalizer, page 71 

ARF Analog Rank Filter, page 81 

AARF Adaptive Analog Rank Filter, page 84 

ARS Analog Rank Selector, page 88 

AQDEF Analog Quantile Density Filter, page 95 

15 AQDOF Analog Quantile Domain Filter, page 95 

AQVF Analog Quantile Volume Filter, page 95 
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Selected Notations and Where They First Appear 
^^D) 

input-output characteristic of a discriminator, Eq. (3) 

Heaviside unit step function, Eq. (4) 

^^^^ Dirac delta function, Eq. (4) 

5 OdJ^ \d{D) .... input-output characteristic of a differential discriminator (probe), Eq, (6) 

j^^""^ volume integral as defined by Eq. (8) 

C= -^ir(x;a,i) cumulative distribution function, Eq. (9) 

ry fK{^;B.,t) density function, Eq. (10) 

J ^<?(a;^) gth quantile for Fi^(x;a,i), Eq. (15) 

ry 10 ga{z;a,/3) convolution transform of the density function f^iz; a), Eq. (19) 

fqi^^i) quantile density, Eq. (24) 

g Rg{a.,t) quantile volume, Eq. (25) 

jj: ^ ■ {■■■)t time averages on a time scale T, Eqs. (32) and (33) 

D {M:^K}^{D,t) weighted mean of K with respect to x, Eq. (41) 

■^^ ('"^R spatial averaging with the test function /r(x), Eq. (47) 

/r [D - x(t)] macroscopic threshold density, Eq. (46) 

^K{'D,t) modulated threshold density, Eq. (52) 

{^xK}r(D,t) mean at reference threshold, Eq. (53) 

^(^^'^) amplitude density, Eq. (54) 

^(D'^) counting density, Eq. (55) 

'^(^^^) counting rates, Eq. (56) 

RC'm time impulse response function, Eq. (57) 

0 estimator of differences in quantile domain between the mean at reference 

threshold and the time average, Eq. (63) 

25 Cic(D,i) modulated cumulative threshold distribution, Eq. (64) 

CK,r [x(t), t] signal x(i), rank normalized with respect to the reference distribution 

CKA'D^t), Eq. (81) 

Knm^ {Kr-)^ {{Kr"^)^y' , Eq. (84) 
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Qab{t; q) . . . .estimator of differences between distributions Ca{D, t) and Cb{D, t), Eq. (92) 

statistic of a type of Eqs. (95) and (97) for comparison of two distributions 

Pg{t) probability that a value drawn from the first sequence is q times larger than the one 

drawn from the second sequence, Eq. (100) 

b(P]t,n[^)) threshold averaged instantaneous density for a continuous ensemble of 

variables, Eq. (120) 

-B(D; t, n{fi)) . . threshold averaged instantaneous cumulative distribution for a continuous 

ensemble of variables, Eq. (121) 

CK(D;i,n(/i)) modulated threshold density for a continuous ensemble, Eq. (122) 

Ci<(D;i, n(^)) . . .modulated cumulative distribution for a continuous ensemble, Eq. (123) 

CK{D;a.,t) modulated threshold density for a scalar field, Eq. (131) 

Ci^(D;a,i) modulated threshold density for a vector field, Eq. (138) 

CkCD; a, t, n{ii)) . . modulated threshold density for an ensemble of vector fields, Eq. (139) 
{MxK}r,A(D; a, t) mean at reference threshold for a vector field input variable, Eq. (140) 
Sq(J^; a, t) quantile domain factor, Eq. (144) 
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DETAILED DESCRIPTION OF THE INVENTION 

The Detailed Description of the Invention is organized as follows. 

In Section 1 (p. 34) we identify the general form of a variable which is subject to analysis 
by this invention, and provide several examples demonstrating the convenience of such a 
5 general representation. 

In Section 2 (p. 35) we use the example of a single scalar variable to describe the basic 
elements of the analysis system adopted in this disclosure, and introduce the discriminators 
and probes as the sensors of such a system. The example of a single scalar variable is 
used to illustrate that the use of discriminators and probes enables us to reformulate many 
10 algebraic problems of the conventional analysis of variables as geometrical problems in the 
threshold space. The particular continuous fashion in which this geometrical extension is 
performed enables the solution of these problems by methods of differential geometry. 

In Section 3 (p. 39) we describe some of the exemplary discriminators and the respective 
probes. 

15 In Section 4 (p. 40) we introduce the normalized scalar fields in the meaning adopted 
for the purpose of this disclosure as the density and cumulative distribution functions 
in the threshold space. We provide a tangible example of how the usefulness of these 
objects extends beyond making details of the analysis intuitive and more available for 
human perception by analogy with the ordinary space. 

20 In Section 5 (p. 42) we provide several examples of equations which reflect the geomet- 
rical properties of the threshold distributions, and are later used for development of various 
practical embodiments of AVATAR. In particular, the definitions of the quantile density, 
domain, and volume are given along with the explanatory examples. 

Section 6 (p. 47) contains a brief additional discussion of possible relationships between 

25 the input and the reference variables. 

In Section 7 (p. 49) we give an introduction to a more general definition of the modulated 
threshold densities by analyzing an example of the threshold crossing density, a quantity 
which cannot be defined for digitized data. 
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In Section 8 (p. 51) we generalize the result of Section 7 by introducing the modu- 
lated threshold densities and the weighted means at thresholds. Along with explanatory 
examples, we show that the weighted mean at reference threshold is indeed a measurement 
of the input variable in terms of the reference variable. We also outline an approach to 

5 computation of the mean at reference threshold by analog machines. 

In Section 9 (p. 55) we interpret the process of measurement by real analog machines 
as a transition from the microscopic to macroscopic threshold densities, or as threshold 
averaging by a probe. Thus we introduce the main practical embodiments of AVATAR as 
the modulated threshold density (Eq. (52)) and the mean at reference threshold (Eq. (53)), 

10 along with the specific embodiments of the amplitude (Eq. (54)) and counting (Eq. (55)) 
densities, and the counting rates (Eq. (56)). We also provide a simplified diagram of 
a continuous action machine implementing the transformation of the multivariate input 
variable(s) into the modulated threshold densities. 

In Section 10 (p. 58) we consider a specific type of weighting function, which is a 

15 convenient choice for time averaging in various embodiments of AVATAR. 

In Section 11 (p. 59) we focus on some of the applications of AVATAR for enhancement 
of analysis through geometric interpretation of the results. We give several examples of 
displaying the modulated threshold densities and provide illustrative interpretation of the 
observed results. Among various examples of this section, there are examples of display- 

20 ing the time evolution of the quantile density, domain, and volume. In this section, we 
introduce such practical embodiments of AVATAR as the phase space amplitude density 
(Eq. (60)), the phase space counting density (Eq. (61)), and the phase space counting rates 
(Eq. (62)). We also provide the illustrative examples of displaying these densities and the 
rates. In Subsection 11.1 we give some illustrative examples of continuous action machines 

25 for displaying the modulated threshold densities and their time evolution. 

In Section 12 (p. 65) we provide a practical embodiment (Eq. (63)) of an estimator of 
differences in the quantile domain between the mean at reference threshold and the time 
average of a variable. 
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In Section 13 (p. 66) we provide a practical embodiment of the modulated cumulative 
distribution (Eq. (64)) and describe how the transition from the densities to the cumulative 
distribution functions in various practical embodiments is formally done by replacing the 
probes by their respective discriminators. Even though the multivariate cumulative distri- 

5 bution function is often disregarded as a useful tool for either graphical or data analytical 
purposes (Scott, 1992, page 35, for example), it is an important integral component of 
AVATAR and is used in its various embodiments. 

In Section 14 (p. 67) we develop simple unimodal approximations for an ideal den- 
sity function, that is, the density function resulting from the measurements by an ideal 

10 probe. Although these approximations are of limited usage by themselves, they provide a 
convenient choice of approximations for the rank normalization. 

In Section 15 (p. 71) we introduce several practical embodiments of rank normaliza- 
tion, such as the general formula for the rank normalization with respect to the reference 
distribution CK,r(D, t) (Eq. (86)), normalization by a discriminator with an arbitrary input- 

15 output response (Eq. (88)), and normahzation of a scalar variable by a discriminator with 
an arbitrary input-output response (Eq. (89)). 

In Section 16 (p. 74) we discuss the usage of the rank normalization for comparison of 
variables and for detection and quantification of changes in variables. We provide several 
simplified examples of such usage and describe a practical embodiment of a simple esti- 

20 mator of differences between two distributions (Eq. (92)). In Subsection 16.1 we provide 
additional exemplary practical embodiments of the estimators of differences between two 
time dependent distributions, Eqs. (95) and (97). In Subsection 16.2 we provide an exam- 
ple of the usage of these estimators for comparing phase space densities and for addressing 
an exemplary speech recognition problem. We also give an outline of an approach to im- 

25 plementation of such comparison in an analog device. In Subsection 16.3 we provide an 
embodiment for a time dependent probabilistic comparison of the amplitudes of two signals 
(Eq. (102)). 

In Section 17 (p. 81) we discuss the usage of AVATAR for analog implementation of 
rank filtering. 
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In Section 18 (p. 82) we discuss the two principal approaches to analog rank filtering of 
a single scalar variable enabled by AVATAR: (1) an explicit expression for the output of a 
rank filter (Subsection 18.1), and (2) a differential equation for the output of such a filter 
(Subsection 18.2). In Subsection 18.1, we also describe a practical embodiment (Eq. (105)) 
for the explicit analog rank filter. 

In Section 19 (p. 84) we briefly discuss the usage of a particular choice of a time weighting 
function in analog rank filters. 

In Section 20 (p. 84) we describe the main embodiment (Eq. (113)) of an adaptive analog 
rank filter. In Subsection 20.1, we also provide an alternative embodiment (Eq. (117)) of 
this filter. 

In Section 21 (p. 87) we extend the definitions of the modulated threshold densities and 
cumulative distributions to include ensembles of variables. We provide the expressions for 
the threshold averaged instantaneous density and cumulative distribution of a continuous 
ensemble, Eqs. (120) and (121), and for the modulated density and cumulative distribution 
of a continuous ensemble of variables, Eqs. (122) and (123). 

In Section 22 (p. 88) we introduce the analog rank selectors, and provide the equations 
for the analog rank selectors for continuous (Eq. (126)) and discrete (Eq. (129)) ensembles. 

In Section 23 (p. 90) we describe the embodiment of an adaptive analog rank filter for 
an ensemble of variables, Eq. (130). 

In Section 24 (p. 90) we introduce the modulated threshold densities for scalar fields, 
Eq. (131). 

In Section 25 (p. 91) we describe the analog rank selectors and analog rank filters for 
scalar fields, Eqs. (133), (134), and (135). In Subsection 25.1, we provide an example of 
filtering monochrome images using a simple numerical algorithm (Eq. (136)), implementing 
an analog rank selector for scalar fields. 

In Section 26 (p. 93) we complete the description of the primary embodiment of 
the AVATAR by generalizing the modulated threshold densities to include vector fields, 
Eq. (138), and ensembles of vector fields, Eq. (139). 

In Section 27 (p. 94)) we provide the description of the mean at reference threshold for 
a vector field input variable, Eq. (140). 
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In Section 28 (p. 95) we describe such important embodiments of AVATAR as the 
analog filters for the quantile density, domain, and volume. These quantities are defined in 
AVATAR for multivariate densities, and thus they are equally applicable to the description 
of the scalar variables and fields as well as to the ensembles of vector fields. 

In Section 29 (p. 97) we provide several additional examples of performance of analog 
rank filters and selectors. 

In Section 30 (p. 102) we provide a summary of some of the main transformations of 
variables employed in this disclosure. 
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1 VARIABLES 



In order to simplify dividing the problem of measurement and analysis of different variables 
into specific practical problems, let us assume that a variable x can be presented as an 
ensemble of vector fields, that is, it can be written as 



where d/x is the weight of the /ith component of the ensemble such that /^d/^n(^) = 
A'', a is the spatial coordinate, and t is the time coordinate. 

Convenience of the general representation of a variable by Eq. (1) can be illustrated by 
the following simple examples. As the first example, consider a vehicle in a traffic control 

10 problem. The position of this vehicle can be described by a single vector variable x — x(t). 
If, however, we are measuring the weight of the vehicle's cargo, it might be more convenient 
to describe it as a scalar field x = x{a.,t), since the weight a; at a given time depends on 
the position vector a. If we now extend our interest to the total weight of the cargo carried 
by different vehicles, this weight is conveniently expressed by an ensemble of scalar fields 

15 X = /^d/in(^) a;^(a, t), where a is the position vector, n{ij) = n(/i) Eili S{iu-i), Xi is the 
cargo capacity (volume of the cargo space) of the i th vehicle, and n{i) is the density (weight 
per unit volume) of the cargo in the i th vehicle. Notice that all the different variables in 
this example can be written in the general form of Eq. (1). 

As a second example, a monochrome image at a given time is determined by the intensity 

20 of the color at location a, and thus it is conveniently described by a scalar field x = a;(a, t). 
A truecolor image is then naturally expressed by a vector field x = x(a, i), where the color 
is described by its coordinates in the three-dimensional color space (red, green, and blue), 
at the position a. We can also consider a "compound" image as a finite or infinite set 
of such images, weighted by the weights n{n). For example, such a compound image can 

25 be thought of as a statistical average of the video recordings, taken by several different 
cameras. 

Additional particular examples of analysis of the variables satisfying the general form 
of Eq. (1) will be provided later in the disclosure. 



5 




(1) 
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2 BASIC ELEMENTS OF SYSTEM FOR ANALYSIS OF VARIABLES 

A system for analysis of variables adopted in this disclosure comprises such basic elements as 
a Threshold Filter, which can be either a Discriminator or a Probe, and an Averaging Filter 
operable to perform either time or spatial average or both time and spatial average. This 

5 system can also include optional modulation and normalization by a Modulating Variable, 
A simplified schematic of such a basic system for analysis of variables is shown in Fig. la. 
This system is operable to transform an input variable into an output variable having 
mathematical properties of a scalar field of the Displacement Variable. The Threshold 
Filter (a Discriminator or a Probe) is applied to a difference of the Displacement Variable 

10 and the input variable, producing the first scalar field of the Displacement Variable. This 
first scalar field is then filtered with a first Averaging Filter, producing the second scalar 
field of the Displacement Variable. Without optional modulation, this second scalar field is 
also the output variable of the system, and has a physical meaning of either an Amplitude 
Density (when the Threshold Filter is a Probe), or a Cumulative Amplitude Distribution 

15 (when the Threshold Filter is a Discriminator) of the input variable. 

A Modulating Variable can be used to modify the system as follows. First, the output 
of the Threshold Filter (that is, the first scalar field) can be multiphed (modulated) by 
the Modulating Variable, and thus the first Averaging Filter is apphed to the resulting 
modulated first scalar field. For example, when the Threshold Filter is a Probe and the 

20 Modulating Variable is a norm of the first time derivative of the input variable, the output 
variable has an interpretation of a Counting (or Threshold Crossing) Rate. The Modulat- 
ing Variable can also be filtered with a second Averaging Filter having the same impulse 
response as the first Averaging Filter, and the output of the first Averaging Filter (that is, 
the second scalar field) can be divided (normalized) by the filtered Modulating Variable. 

25 As will be discussed further in the disclosure, the resulting output variable will then have a 
physical interpretation of either a Modulated Threshold Density (when the Threshold Filter 
is a Probe), or a Modulated Cumulative Threshold Distribution (when the Threshold Filter 
is a Discriminator). For example, when the Threshold Filter is a Probe and the Modulating 
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Variable is a norm of the first time derivative of the input variable, the output variable will 
have an interpretation of a Counting (or Threshold Grossing) Density. 

Let us now describe the basic elements of the analysis system adopted in this disclosure 
in more details, using the measurement a single scalar variable x = x{t) as an example. 

5 We assume the data acquisition and analysis system which comprises the elements 

schematically shown in Fig. lb. Let the input signal be a scalar function of time x{t). This 
input signal (Panel I) is transformed by the discriminator (Panel Ila) into a function of two 
variables, the time t and the displacement D. The latter will also be called threshold in the 
subsequent mathematical treatment. The result of such transformation of the input signal 

10 by the discriminator is illustrated in Panel Ilia. We will provide illustrative examples of 
such a measuring system below. 

The input-output characteristic of the discriminator is described by the continuous 
monotonic function J^az?(^) and is illustrated in Panel Ila. We shall agree, without loss of 
generality, that 

Jkxi^TADix) = 0 

:Fad{0) = 1 ■ (2) 

Although the convention of Eq. (2) is not necessary, it is convenient for the subsequent 
mathematical treatment. One skilled in the art will now recognize that a discriminator 
can thus be interpreted as a threshold integrating filter. The dependence of J^ad on the 
width parameter AD then can be chosen in such way that J^ad approaches the Heaviside 
20 unit step function (Arfken, 1985, p. 490, for example) as AD approaches a suitable limit, 
namely 

lim J^ad{D -x)= 6{D - x) , (3) 

where d{x) is defined as 

I 0 for ;r < 0 

d{x) = £dsS{s) = \ I for x^O , (4) 
I 1 for a; > 0 

25 and d{x) is the Dirac 5-function. When the functional form of the discriminator is the 
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Heaviside unit step function 9{x], such discriminator will be called an ideal discriminator. 
Some other exemplary functional choices for discriminators will be discussed further. 

As an illustration, consider the following example. Imagine that the signal x{t) in 
Panel I is an electrical current which we measure with the ammeter of Panel Ila. The 
scale of the ammeter is calibrated in units of current, and D is our reading of this scale. 
Then J^ad(-D - x) can be interpreted as the probability that our measurement (reading) 
D exceeds the "true" value of the input current x, and AD is indicative of the precision 
of the instrument (ammeter). Thus J^^c [D - x(t)] (Panel Ilia) is just such probabiUty 
with respect to a time-varying input signal, and this probability is now a function of both 
threshold and time. Notice that this function 

z^J'AD[D-x{t)]=f(t,D) (5) 
is represented in a three-dimensional rectangular coordinate system by a surface, which 
is the geometric locus of the points whose coordinates t, D, and z satisfy Eq. (5). Thus 
the methodological purpose of the measurements by a means of discriminators and probes 
can be phrased as "raising" the "flat" problem of the analysis of the curve on a plane 
into the "embossed" surface problem of a three-dimensional space. This alone enables new 
methods of analysis of the signal x{t), and allows more effective solutions of the existing 
problems of the prior methods. As a simple analogy, consider the problem of constructing 
four equilateral triangles out of six wooden matches. This task cannot be achieved on a 
plane, but can be easily accomphshed by constructing a tetrahedron in a three-dimensional 
space. 

The output of the discriminator can be differentiated with respect to the displacement 
(threshold) . The same can be achieved by a means of transforming the input signal by a 
differential discriminator, or probe, as illustrated in Panels lib and Illb. The input-output 
characteristic of a probe is coupled with the one of the discriminator by the relation 

dD^F^niD - x) = ~:Fad{D - x) , (6) 
where dc denotes differentiation with respect to the threshold D. As follows from the 
description of a discriminator, it is convenient, though not necessary, to imagine Od^ad to 
be nonnegative. It is simplest to assume that doJ^ADix) has a single extremum at a; = 0, 
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and vanishes at x = ±00. 

As another example, consider a voltage x{t) (Panel I) applied to the vertical deflecting 
plates of an oscilloscope (Panel lib). If the vertical scale of the graticule is calibrated 
in units of voltage, we can imagine dD^AuiD — x) to be the brightness (in the vertical 

5 direction) of the horizontal line displayed for the constant voltage x, with AD indicative of 
the width of this line. Then OdJ^ad [D - x{t)] (Panel Illb) describes the vertical profile of 
the brightness of the displayed line for the time- varying input signal x{t). 

As will be discussed subsequently, the input-output characteristic of a probe can be 
called the threshold impulse response function of the detection system. One skilled in the 

10 art will recognize that a probe can thus be interpreted as a threshold smoothing filter. The 
functional form of the probe will also be called the (threshold) test function in the sub- 
sequent mathematical treatment. Clearly, as follows from Eqs. (4) and (6), the threshold 
impulse response of an ideal detection system, that is, a system employing ideal discrimi- 
nators, is described by the Dirac J-function. Some other exemplary functional choices for 

15 probes will be discussed further. 

For the purpose of this disclosure, we will further refer to the measuring system com- 
prising the non-ideal discriminators and probes as a "real" measuring system. The output 
of such a system, and thus the starting point of the subsequent analysis, is no longer a 
line in the time- threshold plane (as in the case of an ideal system), but a continuous sur- 

20 face (see Eq. (5), for example). Based on the described properties of the discriminators 
and probes, one skilled in the art will now recognize that discriminators and differential 
discriminators effectively transform the input signal into objects with mathematical prop- 
erties of cumulative distribution and density functions, respectively. The main purpose of 
such transformation is to enable differentiation with respect to displacement (threshold), 

25 while preserving, if originally present, differentiability with respect to space and time. If 
the original input signal is not time-differentiable (e.g., the input signal is time-sampled), 
differentiability with respect to time can always be enabled by introducing time averaging 
into the acquisition system. Likewise, differentiability with respect to spatial coordinates 
can be enabled by spatial averaging. 

30 Particular practical embodiments of the discriminators and probes will depend on the 
physical nature of the analyzed signal(s). For example, for a scalar signal of electrical 
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nature, the discriminator can be viewed as a nonlinear amplifier, and the threshold as a 
displacement voltage (or current, or charge). If the incoming signal describes the inten- 
sity of light, then the displacement can be a spatial coordinate z, and the discriminator 
can be an optical filter with the transparency profile described by ^a/?(^)- The differen- 

5 tial discriminator (probe) can then be implemented through the techniques of modulation 
spectroscopy (see Cardona, 1969, for example, for a comprehensive discussion of modula- 
tion spectroscopy). As an additional example, consider the modification of the previously 
discussed current measurement as follows. Imagine that a gray optical filter is attached 
to the needle of the ammeter, and the white scale is observed through this filter. Assume 

10 that when no current flows through the ammeter, the blackness observed at the position 
D on the scale is J^ad{D), with zero corresponding to the maximum intensity of the white 
color, and "1" representing the maximum blackness. Then J^^ad [D - x{t)] will describe the 
observed blackness of the scale for the time- varying signal x{t) (see Panel Ilia in Fig. lb). 
If the profile of the filter were changed into do^ADiD), then the observed darkness of the 

15 scale will correspond to the output of a probe rather than a discriminator (see Panel Illb in 
Fig. lb). Since the mathematical description of any practical embodiment will vary httle, 
if at all, with the physical nature of the measuring system and analyzed signal, we will fur- 
ther use the mathematical language without references to any specific physical (hardware) 
implementation of the invention. It is also understood that all the formulae manifestations 

20 of the embodiments immediately allow software implementation. 



3 EXEMPLARY DISCRIMINATORS AND PROBES 

Given the input x, the value of 6{D~x) is interpreted as the output of an ideal discriminator 
set at threshold D (see Nikitin et al., 1998, and Nikitin, 1998, for example). Thus the value 
of 5{D — x) = d{D — x) is the output of an ideal probe. 
25 Input-output characteristics of some exemplary discriminators and the respective probes 
are shown in Fig. 2. Notice that although we show only symmetric discriminators, asymmet- 
ric ones can be successfully applied for particular tasks, as well as any linear combination 
of the discriminators. A particular mathematical expression describing the input-output 
characteristic of a discriminator is important only for mathematical computations and/or 
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computer emulations. Any physical device can serve as a discriminator, as long as its input- 
output characteristic and the characteristic of the respective probe satisfy the requirements 
for a test function, e.g., Eqs. (2), (3), and (6). 

For those shown in Fig. 2, the mathematical expressions are as follows: 

Gaussian: dj^J^^^^D) ^ ^^e-i^oY, :Fad{D) 

Cauchy: doJ^ADiD) = ^ [l + (^)'] , Tad{D) 

Laplace: dnJ^ADiD) = ^ e"^ , :Fad{D) 

Hyperbohc: dvJ^AD (^) = ^ (e^ + sir) (D) 

4 NORMALIZED SCALAR FIELDS 



= I erfc (^) ; 

4 + i arctan (^) ; 

4ke^-e^ 



= ^1 + 



tanh (^)] . 



(7) 



Since the invention employs transformation of discrete or continuous variables into objects 
with mathematical properties of space and time dependent density or cumulative distribu- 
tion functions, these properties and their consequent utilization need to be briefly discussed. 
We will further also use the collective term normalized scalar fields to denote the density 
and cumulative distribution functions. Note that the term "time" is used as a designation 
for any monotonic variable, continuous or discrete, common to all other analyzed vari- 
ables, which can be used for sequential ordering of the measurements. Thus "space" is all 
the remaining coordinates which are employed to govern the values of the input variables. 
The term "threshold space" will be used for the coordinates describing the values of the 
variables. For the purpose of this disclosure only, we will further also use the term "phase 
space" , which will be understood in a very narrow meaning as the threshold space employed 
for measuring the values of the variable together with the values of the first time derivative 
of this variable. 

Let us further use the notation for a volume integral as follows: 

y"^d"r /(r) = £jir, . . . J^Jr^ /(r) , (8) 
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where x = (xi, . . . , a;^) and r = (ri,...,r„) are n-dimensional vectors. This definition 
imphes cartesian coordinates, which we also assume in further presentation. If the subse- 
quent equations need to be re- written in curvilinear coordinates (e.g., for the purpose of 
separation of variables) , this can be done by the standard transformation techniques. Refer 
to Arfken, 1985, Margenau and Murphy, 1956, or Morse and Feshbach, 1953, for example, 
for a detailed discussion of such techniques. 

Now let F7i-(x;a,t) be a space and time dependent cumulative distribution function, i.e., 



Phase space density (see Nicholson, 1983, for example) in plasma physics and probability 
density (see Davydov, 1988, and Sakurai, 1985, for example) in wave mechanics are text- 
book examples of time dependent density functions. Another common example would be 
the spectral density acquired by a spectrometer with spatial and temporal resolution (see 
Zaidel' et al., 1976, for example). In Eqs. (9) and (10), the subscript K denotes functional 
dependence of Fk and on some space and time dependent quantity (variable) K. That 
is, although these equations hold for any given space and time, the shape of fK (and, as a 
result, the shape of Fk) might depend on K. The particular way in which such dependence 
is introduced and utilized in this invention will be discussed further in the disclosure. 

The usefulness of density and cumulative distribution functions for analysis of variables 
extends beyond the fact that the geometric representation makes details of the analysis 
intuitive and more available for human perception by analogy with the ordinary space. It 
also lies in one's ability to set up equations better describing the behavior of the variables 
than the algebraic equations of the prior art. As has been discussed previously, for example, 
a level line of the cumulative distribution function of a scalar variable in the time-threshold 
plane corresponds to the output of an order statistic filter. It is an easy-to-envision simple 
geometric image, having many analogies in our everyday experience (e.g., topographical 
maps showing elevations). Since it is a curve of the plane, it is completely determined 




where /jf(x;a, i) is a density function, i.e., /^(x;a, i) > 0, and 




(10) 
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by an algebraic equation with two variables, F{x, y) = 0, with the easiest transition from 
the implicit to the explicit form as a well-known differential equation (see Bronshtein and 
Semendiaev, 1986, p. 405, Eq. (4.50), for example): 

y'{x) = - ^) . (11) 

Fl^{x,y{x)) ^ ^ 

As another example, imagine that the function /;f(x;a,t) in a traffic control problem 
describes the density of cars (number of cars per unit length of the road). Then the 
properties of this density function might be analogous to the properties of the density of 
fluid. They will, for example, satisfy the continuity equation (Arfken, 1985, p. 40, for 
example). Then fluid equations will be the most appropriate for the description of the 
properties of the density function of the traffic problem. Specific applications of the density 
functions will, of course, depend on the physical nature (or applicable physical analogy) of 
the variables involved, that is, on the nature of t, a, and x, and on the constraints imposed 
on these variables. In the next subsection, we will provide several examples of the equations 
involving the threshold distribution and density functions, which might be of general usage 
for analysis of variables. 



5 RANK NORMALIZATION, RANK FINDING, 
AND RANK FILTERING 

Using the above definitions for the normalized scalar fields, we can now introduce several 
examples of additional transformations of variables. Some of these equations will be used 
further in the disclosure for development of various practical embodiments of AVATAR. Let 
us first identify rank normalization and filtering of variables as follows. 

Let us consider a new (dimensionless scalar) variable (field) y{a., t) defined as 

y(a,t) = F;,[x(a,i);a,i]= |_'2''Vr/K(r; a, t) , (12) 
where x(a, t) is an arbitrary variable. Apparently, 0 < y(a, t) <1. Since Ff^ (x; a, t) > 0, 
where denotes a partial derivative with respect to Xi, Eq. (12) defines rank normalization 
of the variable x(a, t) with respect to the reference distribution Fk- Rank normalization 
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transforms a variable into a scalar variable (or scalar field, if the transformed variable is a 
field variable), the magnitude of which at any given time equals to the value of the reference 
distribution evaluated at the value of the input variable at a this time. Thus Fk provides 
a (nonlinear) scale for measurement of x(a,£). We will discuss rank normalization in more 
5 details later in the disclosure. 

Differentiating Eq. (12) with respect to time leads to yet another defining equation for 
an analog rank normalizer (ARN) as follows: 

y = dt Fk{x] a, t) + {x ■ Vx) Fk{x; a, t) , (13) 

where dt denotes a partial derivative with respect to t, and x • Vx in cartesian coordinates 
10 is simply 

x■V^ = J2^^^^^■ (14) 



Let us now introduce another transformation of variables, which can be interpreted as 
rank filtering. By definition, is the qth quantile of F^(x;a,£) when 

FK^Xg] a.,t) = J ' d^r /^^ (r; a,t) = q = constant , (15) 
where 0 < § < 1 is the quantile value. Note that Eq. (15) describes a simple surface in the 
threshold space. When the variable x is a scalar, that is, x{a,t) = rc(a,i), this surface is a 
point on the threshold fine, and thus (as has been previously discussed) (a;,; a, t) = g de- 
scribes a level line in the time-threshold plane. Taking the full time derivative of FKix^; a, t) 
allows us to rewrite Eq. (15) in differential form for a family of equiquantile surfaces in the 
threshold space as 

dt Fk{x,- a, t) + (x, • Vx) FK{xq- a, f) = 0 . (16) 

We can further introduce some constraints on x^, such as constraints on the direction of 

x^, or other subsidiary conditions. For example, if we allow only the nth component of Xq 

to depend on time, the time derivative of this component can be written as 

_ dtFK{xg-a,t) 
^'"""^..„i^i.(x,;a,£)' (17) 

assuming that all derivatives in this equation exist. Eq. (17) thus defines an analog rank 
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filter (ARF). The latter can also be called analog order statistic filter (AOSF), or analog 
quantile filter (AQF). Note that even though we adopt the existing digital signal process- 
ing terminology such as "order statistic filter" , this is done by a simple analogy between 
the geometric extension of the AVATAR and the definitions in the discrete domain. Our 
5 definitions cannot be derived from the algebraic equations defining order statistic filtering 
in digital signal processing. Note also that a particular form of Eq. (17) depends on the 
nature of constraints imposed on x^. For the important special case of the input variable 
as an m-dimensional surface, that is, a scalar field x = z{a.,t), where a = (ai,...,a^), 
Eq. (17) reads as 

In numerical computations, Eqs. (17) and (18) can be considered to be modifications of 
the Newton- Raphson method of root finding (see Press et al., 1992, for example, and the 
references therein for the discussion of the Newton- Raphson method). 

When the distribution function does not depend on time explicitly, Fk = Fk{z\ a), we 
15 can introduce an explicit parametric dependence of the density function on some /?, for 
example, through the convolution transform 

g^{z- a, (5) = l^de MP - e) Mz; a) , (19) 

such that ga{z;a,i3) approaches fK{z]a.) as a approaches a suitable limit, for example, 
when 0a(,/3 — e) approaches the Dirac 5-function S{/3 — e). Then the equality 

20 /_'^^''^^d£5a(^;a,/3) = g (20) 

leads to the equation for an analog rank finder, or analog rank selector, as follows: 

--fa 31- J't-''<i,d,g^{,;^,P) 
d/s"''"'"^- <,.h(a,/?);a,/3] ' <'l) 

Introducing parametric dependence through the convolution transform will later be shown 
to be convenient in rank selectors for an ensemble of time dependent variables, since 
25 the parameter p can be chosen to be the time itself. Clearly, there are plenty of alter- 
natives for introduction of such parametric dependence. For example, one can choose 
/K(e;a, a) = 0(a) //<:(£:; a) with (f){a) such that lima^oo <f>{oi) = 1. As an illustration, the 
choice (/){a) = 1 — (1 — g)e~" leads to 
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da ^ [l-(l-9)e-"]/K(^5;a) ' ^''^ 
where 2q(a, a) rapidly converges to the "true" value of ^^^(a). Notice again that the def- 
initions of the analog rank selectors require the existence of fx, that is, the threshold 
continuity of F^, and thus cannot be introduced in the digital domain. 

Finally, the threshold continuity of the distribution function allows us to write, for 
both time dependent and time independent scalar fields, an explicit expression for the gth 
quantile of Fk{x] a, t) as 

a;,(a, t) = j'^dr r fK(r; a, t) 6 [i^j^ (r; a, t) - q] . (23) 

The derivation and properties of this equation for rank filtering will be discussed later in 
the disclosure. Further, we will also provide a means of evaluating this expression by analog 
machines. 

Let us also define another type of rank filtering, which, unlike the rank filtering defined 
earlier, is applicable to multivariate variables and does not (and cannot) have a digital 
counterpart, as follows: 

£^d"r /^(r; a, t) d [/^(r; a, t) - /^(a, t)] = q = constant , (24) 

where fg{a,t) is quantile density. Since the density function /;<:(x;a, t) vanishes at Xn = 
±oo, the surface in the threshold space defined by Eq. (24) encloses the series of volumes 
(regions in the threshold space) such that /i<r(x; a, t) > /g(a, t), and the integral of /i<'(x; a, t) 
over these volumes is equal to q. We shall call this series of regions in the threshold space 
the quantile domain. Notice that the left-hand side of Eq. (24) is non-increasing function 
of fq{3i,t) > 0. Later in the disclosure, we will provide a means of finding fq{a.,t) by 
continuous action machines. 

We shall designate the total volume enclosed by the surface defined by Eq. (24) as 
quantile volume i?^, which can be computed as 



,(a,i) = ^Vr^[/;,(r;a,t)-/,(a,i)] . (25) 

Notice that the quantile density indicates the value of the density likely to be exceeded, and 
the quantile volume gives the total volume of the highest density regions in the threshold 
space. As an example, consider the density of the cars in a city. The median density 
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will indicate the degree of congestion on the roads, providing the number of cars per unit 
length of a road (or, inversely, "bumper-to-bumper" distance) such that half of the traffic is 
equally or more dense. Then the median domain will indicate the regions (stop lights and 
intersections) of such congested traffic, and the median volume will give the total length 
of the congested roads. As another simple example, consider the price (amount per area) 
of the land in some geographic region. The median density will be the price such that the 
total amount paid for the land equally or more expensive will be half of the total cost of the 
land in the region. Then the median domain will map out these more expensive regions, 
and the median volume will give the total area of this expensive land. Notice that even 
though in the latter example the "median density" is price (it has the units of amount per 
area), it is not the "median price" in its normal definition as "the price such that half of 
the area is less (or more) expensive". Later in the disclosure, we will provide a means of 
computation of the quantile domain and volume by continuous action machines. 

Since the distribution and density functions Fk and fx depend on the quantity K, 
comparison of these functions with respective Fk' and Jk' for a different quantity K' will 
provide a means for assessment of K and K' in their relation to the reference variable, 
that is, to the variable for which the distribution and density were computed. For example, 
the equality Fk = Fk' when K K' will indicate that even though K and K' are not 
equal, they are equivalent to each other in their relation to the reference variable, at least 
under the conditions of the conducted measurement. When the quantity K represents 
the reference variable in some way so that the behavior of K reflects some aspects of the 
behavior of the reference variable, the reference variable should be considered a component 
of the measured signal or process rather than a part of the acquisition system. In this case, 
we shall designate the reference variable as the input variable and call the quantity K an 
associated variable, or an associated signal. When such interdependence between K and x 
is not only implied or suspected, but defined in some particular way, we will also call x the 
input variable, and K a property of the input variable. 

The particular way in which the dependence of the density and the cumulative distri- 
bution functions on K are introduced in this invention can be interpreted as measuring 
the input variable K in terms of the rate of change of the reference variable x at a certain 
threshold D. The details of this interpretation will be provided later in the disclosure. In 



46 



the next subsection, we will briefly discuss some aspects of the relationship between the 
input and the associated variables. 

6 RELATIONSHIP BETWEEN INPUT AND ASSOCIATED 
VARIABLES 



In order to implement comparison of variables of different natures, it is important to have 
5 a reference system, common for all variables. Naturally, time (t) is one of the coordinates 
in this reference system, since it equally describes evolution of all measured variables. Time 
can also serve as a proxy for any monotonically ordered index or coordinate. For brevity, 
we shall call the remaining coordinates governing the values of the variables the spatial 
coordinates (a). Time and space are the third and second arguments, respectively, in the 
10 dimensionless object, cumulative distribution function, we are to define. We will call the 
first argument of the cumulative distribution the threshold, or displacement (D), and the 
units of measurements of this argument will be the same as the units of the input (reference) 
variable. 

There is plenty of latitude for a particular choice of an associated variable K, and 
15 a variety of ways to introduce the coupling between K and the input variable. For the 
purpose of this disclosure, it would be of little interest to us to consider an a priori known 
K other than K = constant. Thus K must not be confused with the part of the acquisition 
system and rather should be considered a component of the measured phenomenon. The 
relationship between K and the input variable can be of a deterministic nature, such as 
20 mathematical transformation, or a result of physical dependence. For example, the reference 
variable can be the instantaneous air temperature, and the associated variable can be the 
instantaneous humidity. Then we might be interested in measuring the dependence of the 
humidity on the temperature variations at a given temperature. Or, the reference variable 
can be the total population of the psychiatric wards in the United States, and the associated 
25 variable the Dow Jones Industrial Average. Then one might try to investigate how the 
rate of change in the mental health of the nation affects the economic indicators. One 
skilled in the art will recognize that such dependence between the input and the reference 
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variable is most naturally described in terms of their joint distribution. However, such 
a joint distribution will be a function of the threshold coordinates of both the input and 
the reference variables. Thus a different input variable will require a different threshold 
space for the description of its dependence on the reference variable. In order to enable 
5 comparison between input variables of different natures, we would desire to use the reference 
system common to both input variables, that is, the threshold coordinates of the reference 
variable. 

As we will try to illustrate further, the choice of both reference and associated variables 
should be based on the simplicity of treatment and interpretation of the results. For our 

10 illustrative purposes, pursuant to easy and useful interpretation, we introduce nonconstant 
associated variables only as norms of the first two time derivatives of the input signal, 
|x| and |x|. Our choice of coupling as modulation, as further described in detail, is based 
solely on the immediate availability of physical interpretation of the results in the case 
of iiT = |x|. For example, the cases K = constant and K = |x| relate to each other as 

15 the charge and the absolute current in electric phenomena. This coupling (as modulation) 
allows the input (reference) variable to provide a common unit, or standard, for measuring 
and comparison of variables of different nature. This coupling also enables assessment of 
mutual dependence of numerous variables, and for evaluation of changes in the variables and 
in their dependence with time. For example, dependence of economic indicators on social 

20 indicators, and vice versa, can be analyzed, and the historical changes in this dependence 
can be monitored. Different choices of associated variables, however, may benefit from 
different ways of coupling. 

For the purpose of this disclosure, we assume that continuous is synonymous to dijfer- 
entiable. Whenever necessary, we assume a continuous input variable x(a, t). When the 

25 definition of a particular property requires continuity of derivatives of the input variable, 
such continuity will also be assumed. For instance, the definitions of the densities for the 
threshold accelerations and for the phase space threshold crossing rates of a scalar variable 
x{t) will require continuity of the first time derivative x{t) of the input variable. Of course, 
all the resulting equations are applicable to digital analysis as well, provided that they are 

30 re- written in finite differences. 
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T THRESHOLD DENSITY FOR COUNTING RATES 
OF SINGLE SCALAR VARIABLE 

As an introduction to a more general definition, let us consider a single scalar continuous- 
time variable (signal) x{t), and define a time dependent threshold density for this signal's 
counting (threshold crossing) rates. First, we notice that the total number of counts, i.e., 
5 the total number of crossings of the threshold D by the signal x{t) in the time interval 
0 < t < T, can be written as (see Nikitin et al, 1998, and Nikitin, 1998, for example) 

N{D) = J2l^dt5it-ti), (26) 

where d{x) is the Dirac (5-function, and t, are such that x{ti) = D for all i. Using the 
identity (Rumer and Ryvkin, 1977, p. 543, for example) 

*-/wi = i:^|i(^. (27) 

we can rewrite Eq. (26) as 

N{D)^ j^dt\x{t)\6[D-x{t)] , (28) 
where the dot over x denotes the time derivative. In Eq. (27), \f'{x)\ denotes the absolute 
value of the function derivative with respect to x and the sura goes over all Xi such that 
15 f{x^) = a. Thus the expression 

^{D) = ^ £dt\xit)\5[D~xit)]^ 

^ /-i^ ^^'^T"'^ \^i^)\nD-x{s)] (29) 
defines the counting, or threshold crossing, rate. 

Fig. 3 illustrates the counting process for a continuous signal. This illustration utilizes 
the fact that 5{x) = ^^(a^), where e{x) is the Heaviside unit step function, and thus 
20 9[x{t) -D]= i;{t) 5[x{t) - D] by differentiation chain rule. 

Replacing the rectangular weighting function in Eq. (29) by an arbitrary time window 
h{t), the rate of crossing of the threshold D by the signal x{t) can be written as the 
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convolution integral 

nH{D,t) = l^dshit- s)\x{s)\S[D -x{s)] , (30) 

where the time weighting function h{t) is such that 

J^Jish{s)^l, (31) 

and is interpreted as a moving (or sliding) window. It is simplest, though not necessary, 
to imagine h{t) to be nonnegative. Notice that now the threshold crossing rate (Eq. (30)) 
depends on time explicitly. Note also that in Eq. (30) this rate is measured by means of an 
ideal probe do^AniD - x) =^ 5{D - x) with time impulse response h{t). 

If T is a characteristic time, or duration of h{t), we will use shorthand notations for the 
integral 

J2dsh{t-s-T) ••• = (••.)^ (32) 
= {■■■)t. (33) 
These equations define a time average on a time scale T. We will consequently use the 
notations of Eqs. (32) and (33) for both continuous and discrete time averages. The no- 
tation of Eq. (32) will be used whenever the particular choice of the weighting function is 
important. 

Noticing that j^^(lDnh{D,t) = {\x\)j,, we can now define the counting (threshold 
crossing) density as 

. {\x\6iD -x))j. , ^ 

mr ' ^''^ 

where we used the shorthand notation of Eq. (33). 

The meaning of the above equation can be clarified by its derivation from another 
simple reasoning as follows. Note that a threshold crossing occurs whenever the variable 
has the value D, and its first time derivative has a non-zero value. Then the density of 
such events is expressed in terms of the joint density of the amplitudes of the variable and 
its time derivative as 
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r(D t) = /-°°ood^^l-P^I {5{Di-i)5{D-x))^ 



(^S^^dD J^^dDi\Di\5{Di-x) S{D -x)^ {\x\)t 



(35) 



The significance of the definition of the time dependent counting (threshold crossing) 
density, Eq. (34), stems from the importance of zero- crossings, or, more generally, threshold 
crossings, and zero/threshold crossing rates for many signal processing applications. These 
quantities characterize the rate of change in the analyzed signal, which is one of the most 
important characteristics of a dynamic system. The importance of threshold crossing rates 
can be illustrated by the following simple physical analogy: If x{t) describes the location 
of a unit point charge, then 6 {D - x) is the charge density, and thus \x\ 5{D ~x) is the 
absolute current density at the point D. In the next subsection, we generalize the above 
result and provide its additional interpretation. 

8 MODULATED THRESHOLD DENSITIES AND WEIGHTED 
MEANS AT THRESHOLDS 



In order to generalize the above result, let us first analyze the example shown in Fig. 4. 
Consider intersections of a scalar variable a: (t) in the interval [0,T] with the thresholds [Dj], 
15 where Dj+i — Dj -\- AD. The instances of these crossings are labeled as {ti}, ti+i > ti. 
The thresholds {Dj} and the crossing times {ti} define a grid. We shall name a rectangle 
of this grid with the lower left coordinates {ti,Dj) as a Sij box. We will now identify the 
time interval At^ as ti+i - ti if the box Sy covers the signal (as shown in Fig. 4), and zero 
otherwise. 

20 We can thus define the threshold density, modulated (weighted) by the associated vari- 
able K, or simply Modulated Threshold Density (MTD), as 

^ Y.^UoK{ti) 

CK{Dj,t) = lim -r^ ^ . (36) 

^' ^ AD-,0 AD J2AtijK{ti) ^ ^ 
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Utilizing Eq. (27), we can rewrite Eq. (36) as 

For example, K{t) = \x{t) \ leads to the previously described result, that is, to the counting 
density. For K = constant, Eq. (37) reduces to the amplitude density (Nikitin, 1998), 
namely 

6(D, t) = J^Js hit -s)5[D- X{s)] = (S{D- x))^ . (38) 

Notice that the modulated threshold density also formally reduces to the amplitude density 
whenever 

{K6{D~ x))^ = (K)^ {5 (D - x))^ , (39) 
that is, when K{t) and the pulse train 5[D - x{t)] =: \x{t)\-'^ Y.i 5[t - ti) are uncorrelated . 

To further clarify the physical meaning of MTD, let us first use Eq. (27) and rewrite 
the numerator of Eq. (37) as 

{K5{D-x)% = Y.Kt-U)~^^, (40) 

which reveals that it is just a weighted sum of the ratios Kil\x{ti)l evaluated at the in- 
tersections of x with the threshold D. Noticing that the ratio dl?/|a;(i,)| is equal to the 
time interval the variable x spends between D and D + dD at the ith intersection, we shall 
realize that 

is the time weighted mean of K with respect to x at the threshold D, or simply Mean 
at Reference Threshold (MRT). Using this designation (MRT) can further be justified by 
rewriting the middle term of Eq. (41) as 

{K6iD-x)f^ ^ { [J-oo^Dk Dk 5 [Dk - K)] 5 (D - x))^ 
{6{D-x))!^ ~ {5{D-x))i, 
f°^^^ ^ {5iDK-K)5(D-x))^ 

which demonstrates that MRT is indeed the first moment of the density function 
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/{Dk; D, t). Notice that here the threshold coordinates D of the variable x are the spatial 
coordinates of the density function /. Later in the disclosure, we will generalize this result 
to include the multivariate modulated threshold densities. Then, for example, if the refer- 
ence variable yi{t) is the positions of the particles in an ensemble, and K(t) = x(t) = v(i) 
is their velocities, then {MxK}^(D,i) will be the average velocity of the particles at the 
location D as a function of time t. Using Eq, (41), Eq. (39) can be rewritten as 



and understood as the equality, at any given threshold, between the weighted mean of K 
with respect to x and the simple time average of K. That is, if Eq. (43) holds for any 
threshold, the weighted mean of K with respect to a; is a function of time only. Obviously, 
Eq. (43) always holds for K = constant. 

It is very important to notice that although the modulated threshold density given by 
Eq. (37) implies that K never changes the sign from "plus" to "minus" (or vice versa), 
the weighted mean at threshold of a reference variable defined by Eq. (41) is meaningful 
for an arbitrary K, and thus the comparison of Eq. (43) can always be implemented. For 
example, this comparison is implementable for {K)'^ = 0, when the modulated density does 
not exist. 

Eq. (43) signifies that a simple comparison of the weighted mean at threshold {M^K}!^ 
with the simple time average {K)!^, that is, comparison of a modulated density with un- 
modulated (amplitude), will indicate the presence of correlation between the reference and 
the associated variables as a function of threshold (and time) on a given time scale. In 
other words, the equality of Eq. (43) holds when, at a given threshold D, the values of the 
variable K are uncorrelated with the time intervals the variable x spends between D and 
D + dD. As a simplified example, consider K{t) as a clipped x{t), that is, as 



In this case, Eq. (43) will hold for D < Dq, and will generally fail for D > Dq. 

As another example, let the signal x{t) represent a response of a detector system to a 
train of pulses with high incoming rate, Poisson distributed in time. The high rate might 
cause high order pileup effects in some low energy channels of the detector. Then, as follows 



(43) 




Do for x{t) < Do 
x(t) for x{t) > Do 



(44) 
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from (Nikitin, 1998), the amplitude and the threshold crossing densities for such a signal, 
measured in these channels, will be identical. Thus the channels afflicted by the pileup 
effects can be identified by comparing the counting rate with the amplitude distribution in 
different channels. 

5 As an opposite extreme of Eq. (43), the mean at reference threshold can be a function 

of threshold only, and thus K would be completely determined by the reference variable. 
As an illustration, consider a linear reference variable x = at, a > 0, and thus d{D — x) = 
H (f ~ 0 ' ^^^ich leads to the MRT as 



10 which depends only on the threshold. One skilled in the art will now recognize that a simple 
reincarnation of Eq. (45) in a physical device would be an ideal oscilloscope (that is, a precise 
oscilloscope projecting an infinitesimally small spot on the screen), where the reference 
variable x = at is the voltage across the horizontal deflecting plates, and the MRT is the 
vertical position of the luminescent spot on the screen a.t x = D. Thus the measurement 

15 of the MRT is indeed "measuring the input variable K in terms of the reference variable 
x" . For an arbitrary reference variable, the MRT is the average of these vertical positions 
weighted by the time intervals the reference variable spends at the horizontal position 
D. Since the afterglow of the luminophor coating of the screen conveniently provides the 
(exponentially forgetting) time averaging, the MRT can be measured as the average vertical 

20 deflection, weighted by the brightness of the vertical line at the horizontal deflection D. 
In the next three subsections we will clarify how this idealized example relates to the real 
measurements. 

Fig. 5 provides another example of using modulated densities for measuring the input 
variable K in terms of the reference variable x. Notice that the amplitude densities (center 

25 panels) of the fragments of the signals xi [t) and X2 {t) shown in the left-hand panels of 
the figure are identical. Notice also that the modulating signals Ki{t), K2{2), and Kait) 
are identical for the respective modulated densities of the signals xi{t) and X2{t), while 
the modulated densities are clearly different. Thus even though the amplitude densities 
and the modulating signals are identical, different reference signals still result in different 

30 modulated densities. 
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Although the first argument in the density function ck{D, t) is always a threshold value, 
we will call the modulated densities for K = constant, K = \x\, and ilT = |:r|, for brevity, 
the amplitude, counting (or threshold crossing), and acceleration densities, respectively. We 
now proceed with the general definitions of the multivariate modulated threshold densities 
and the weighted means at thresholds. 

9 MULTIVARIATE THRESHOLD DENSITIES AVERAGED WITH 
RESPECT TO TEST FUNCTION 

Time averaging with a time weighting function h{t) signifies transition from microscopic 
(instantaneous) densities (5(D - x) to macroscopic (time scale T) densities ((5(D - x))^. 
Carrying out the same transition from microscopic to macroscopic threshold domain can be 
done by a means of averaging with respect to a test Junction /r(x), a standard approach in 
such fields as electrodynamics (Jackson, 1975, Section 6.7, for example) or plasma physics 
(Nicholson, 1983, Chapter 3, for example). We thus can define a (multivariate) macroscopic 
threshold density as a threshold average on a scale R, namely as 

{S{B - r))i = |_Vr /r(x - r) 6{B - r) = /r(D - x) , (46) 

where R is a characteristic volume element. We will assume /r(x) to be real, nonzero 
in some neighborhood of x = 0, and normalized to unity over all space. It is simplest, 
though not necessary, to imagine /r(x) to be nonnegative. Such threshold averaging with 
respect to a test function reflects finite amplitude resolution of data acquisition systems. 
For digitally recorded data, the lower limit for the characteristic volume is the element of 
the threshold grid. We will further use the shorthand notation 

/^d"r/K(x-r) ■■■ = {---}i (47) 
to denote the spatial averaging with the test function /r(x). 

For hardware devices, the choice of /r(x) is dictated by the threshold impulse response 
of the discriminators (probes) (Nikitin, 1998, Chapter 7, for example). In software, this 
choice is guided by computational considerations. For the purpose of this disclosure, it is 
convenient to assume that the total threshold impulse response function is the product of 



55 



Fi*rTt Named* Inventor: Nildtin 
US Express- EL388608922US 

the component impulse responses, that is, it can be written as 

M^) = f[dD^T^DX^,), (48) 

1=1 

where capital pi denotes product, as capital sigma indicates a sum, that is, 

flfi=flf2f3 ... fn- (49) 
i=l 

5 

Unless otherwise noted, the subsequent computational examples employ Gaussian test 
function, namely 

/K(x) = n9.,^.,(..) = ^^exp[-|:(^)] . (50) 
Notice that, for the Gaussian test function, 
10 rS'M"^) = = U^AD.ix,) = 2-" nerfc (^) , (51) 

where erfc(a;) is the complementary error function (Abramowitz and Stegun, 1964, for 

example). .?-r(x) should be interpreted as threshold step response. 

Fig. 6 illustrates an optical threshold smoothing filter (probe). This probe consists of a 

point light source 5' and a thin lens with the focal length /. The lens is combined with a 
15 gray optical filter with transparency described by /2/(x). Both the lens and the filter are 

placed in a XOY plane at a distance 2/ from the source S. The lens-filter combination can 

be moved in the XOY plane by the incoming signal r so that the center of the combination 

is located at in this plane. Then the output of the filter is proportional to the intensity 

of the light measured at the location D = {D^,Dy) in the D^-O-Dy plane parallel to the 
20 XOY plane and located at the distance R from the image 5" of the source S (toward the 

source). That is, the output of this filter can be described by /r(D — r). 

When a test function is employed for threshold averaging, Eq. (37) can be rewritten 

for the multivariate modulated threshold densities as Threshold- Time Averaged Density, 

namely as 

- ____ . 
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As for a scalar variable, in this disclosure we shall call the multivariate modulated 
threshold densities for K = constant the amplitude densities, and the densities for 
K = [YJl^i{xJ /\D^YY the (multivariate) counting densities. The amplitude density thus 
indicates for how long the signal occupies an infinitesimal volume in the threshold space, 
and the counting density indicates how often the signal visits this volume. For example, 
the amplitude density will indicate the number of cars at a certain intersection (that is, at 
the intersection positioned at D) at a certain time, and the counting density will describe 
the total traffic flow through this intersection at a given time. Carrying out the threshold 
averaging in Eq. (41), we can write the equation for the mean at reference threshold as 

{ /r [D - x(s)] >^ 

= (K(5) nr=i^A^AA [A-x.(g)])y , 
{ nr=i dv.^AD. [A - x,{s)] )y • ^^'^^ 

Notice that, unlike in the definition of the modulated threshold density, the variable K no 
longer has to be a scalar. In the traffic example above, it is obvious that the same traffic 
flow (number of cars per second) can be achieved either by low density high speed traffic, 
or by high density low speed traffic. The ratio of the counting rates and the amplitude 
density will thus give us the average speed at the intersection. If the variable K is different 
from the speed - for example, it is the rate of carbon monoxide emission by a car - then 
the mean at reference threshold will indicate this emission rate at location D, and it may 
(for example, because of the terrain or speed limit) or may not depend on the location. 

In this disclosure, we assume the validity of Eq. (48), and thus the explicit expressions 
for the amplitude density b(D, t) and the counting density r(D,t) are as follows: 

6(D, t) = (^f[ Od^^Pad^ [A - x,{s)] j (54) 

for the amplitude density, and 

( [Ig]' nr=i ^i^.^aa [a - x,{s)] ) 

r(D, t) = IJT (55) 

for the counting density. In Eq. (55), the numerator in the right-hand side is proportional 
to the counting rate. Explicitly, the expression for the counting rates reads as 
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*) = { E [1^] " n AA dn.T^o, (A - a:,(.s)] }^ . (56) 

A type of problem addressed by Eq. (56) might be as follows: How often does a flying 
object cross a radar beam? 

Fig. 7 shows a simplified diagram illustrating the transformation of an input variable into 

5 a modulated threshold density according to Eq. (52). The sensor (probe) of the acquisition 
system has the input-output characteristic /r^ of a differential discriminator. The width 
of this characteristic is determined (and may be controlled) by the width, or resolution, 
parameter R^. The threshold parameter of the probe D signifies another variable serving 
as the unit, or datum. In Fig. 7, the input variable x^(i) is a scalar or vector, or a component 

10 of an ensemble. For example, a discrete surface (such as an image given by a matrix) can be 
viewed as a discrete ensemble, being scalar for a monochrome image, and a 3D-vector for 
a truecolor image. The output of the probe then can be modulated by the variable K^{t), 
which can be of a different nature than the input variable. For example, Kf,{t) - constant 
will lead to the MTD as an amplitude density, and K^{t) = |x^(t)| will lead to the MTD 

15 as a counting density/rate. Both the modulating variable if^ and its product with the 
output of the probe iT^ /r^ can then be time-averaged by a convolution with the time 
weighting function h{t;T), leading to the averages (i^^/R^(D - x^))^ and The 
result of a division of the latter average by the former will be the modulated threshold 
density CK^i'D,t). Notice that all the steps of this transformation can be implemented by 

20 continuous action machines. 

10 TIME AVERAGING OF MULTIVARIATE THRESHOLD 
DENSITIES BY RCm IMPULSE RESPONSE FUNCTIONS 



Let us consider a specific choice of a time weighting function as follows: 

^nit) = -f^t-e-^eit). (57) 
This is a response of a circuit consisting of one RC differentiator and n RC integrators (all 
25 time constants RC = T) to a unit step of voltage e(t). Thus we shall call such weighting 



58 



Atty. Docket No 152S2 2 

First Named Inventor- Nikitm 
US Express. EL3S8608922US 



function an RCi^ impulse response. Fig. 8 shows the RCin time weighting functions for 
n = 0 (exponential forgetting), n = 1, and n = 2. 
Differentiation of hn{t) leads to 

K{t)=^^ [K-i{t) - K{t)] (58) 

5 for n > 1, and 

ho{t) = ^[5{t)-h,{t)] (59) 

for n = 0. 

An exponential factor in time weighting functions is ubiquitous in nature as well as 
in technical solutions. In particular, i?Ci„ impulse response time averaging functions are 
10 quite common, and easily implementable in software as well as in various devices. Although 
normally the time weighting function does not need to be specified in detail, an exponential 
factor in the time weighting function allows us to utilize the fact that (e^)' = e^. In partic- 
ular, the relations of Eqs. (58) and (59) allow us to simplify various practical embodiments 
of AVATAR. This will become apparent from further disclosure. 



11 SHAPE RECOGNITION AND 
TIME EVOLUTION OF DENSITIES. 
DISPLAYING DENSITIES AND THEIR TIME EVOLUTION 

15 

As has been mentioned earlier in this disclosure, the main purpose of the analysis of 
variables through their continuous density and distribution functions is twofold: (1) to 
facilitate the perception through geometric interpretation of the results, and (2) to enable 
the analytical description by differential methods. Let us first address the former part 

20 of this goal, that is, the visual presentation of the densities and the interpretation of the 
underlying qualities of the variable based on these observations. 

Let us first notice that the amplitude density at a given threshold D is proportional to 
the time the variable spends around this threshold, and thus is proportional to the average 
inverted (absolute) slope of the variable at this threshold. One might say that the ampli- 

25 tude density is a measure of "flatness" of the signal. For example, the amplitude density 
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generally increases with the increase in the number of extrema at the threshold D. The 
counting density is proportional to the number of crossings, or "visits", of this threshold by 
the variable, and the acceleration density is generally proportional to the density of sharp 
turns (such as extrema and inflection points) of the variable at this threshold. Thus the 

5 acceleration {K = \x\), amplitude {K ^ 1), and counting {K = \x\) densities complement 
each other in a manner necessary for selective shape recognition of a signal, as illustrated 
in Figs. 9 a and 9 b. The left columns of the panels in these figures show the fragments of 
three different signals in rectangular windows. The second columns of the panels show the 
amplitude densities, the third columns show the counting densities, and the right columns 

10 show the acceleration densities for these fragments. These figures illustrate that the ac- 
celeration and counting densities generally reveal different features of the signal than do 
the amplitude densities. For the fragment Xi(t) in Fig. 9 a (the upper row of the panels), 
\x{t)\ = constant, and thus the counting and the amplitude densities are identical. For the 
fragment X2{t) in Fig. 9a (the middle row of the panels), \x{t)\ = constant, and thus the 

15 acceleration and the amplitude densities are identical. 

The example in Fig. 10 shows time dependent acceleration densities, threshold crossing 
rates, and amplitude densities computed in a 1-second rectangular moving time window 
for two computer generated non-stationary signals (Panels la and lb). Panels 2a and 2b 
show the acceleration densities. Panels 3a and 3b show the threshold crossing rates, and 

20 Panels 4a and 4b show the amplitude densities. The signals represent sequences of (non- 
linearly) interacting unipolar Poisson-distributed random pulses, recorded by an acquisition 
system with an antialiasing bandpass RC-filter with nominal passbands 0.5-70 Hz at -3 dB 
level. The sequences of the pulses before interaction are identical in both examples, but the 
rules of the interaction of the pulses are slightly different. These differences are reflected 

25 in the shape of the resulting signals, which can in turn be quantified through the displayed 
densities and the rates. 

As an example of the interpretation of the displayed densities, consider the stretch of 
the signals in the interval 45 through 70 seconds. For both signals, the amplitude density 
(Panels 4a and 4b) is highest at the lowest amplitude, and is approximately uniform at 

30 other thresholds. This is likely to indicate that the signals in this time interval consist of 
relatively narrow tall pulses of comparable amplitude, originating from a flat background. 
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The approximate uniformity of the counting rates (Panels 3a and 3b) between the lowest 
and the highest thresholds confirms the absence of the "secondary" extrema, that is, these 
are single pulses. The increased rates in the intervals 50 to 60 seconds and 65 to 70 seconds 
indicate that there are more pulses per unit time in these intervals than in the intervals 45 

5 to 50 and 60 to 65 seconds. An approximate equality of the acceleration densities (Panels 2a 
and 2b) at the highest and lowest thresholds is likely to indicate that these pulses might 
have sharp onsets and sharp "tails" , in order for the "sharpness" of the peaks to be equal 
to the combined "sharpness" of the onsets and the tails. 

Earlier in this disclosure, we introduced a new type of rank filtering, which is applicable 

10 to analysis of scalar as well as multivariate densities. In Eqs. (24) and (25), we introduced 
the quantile density, domain, and volume. Let us now illustrate how these quantities are 
applicable to the analysis of scalar variables. Panel I of Fig. 11 shows the fragment of 
the signal from Panel la of Fig. 10 in the time interval between 7 and 29 seconds. The 
amplitude density of this fragment is plotted in Panels II through IV. In these panels, the 

15 quartile densities /1/4 (Panel II), /1/2 (Panel III), and /3/4 (Panel IV) are shown by the 
horizontal lines. These lines intersect the density in such a way that the shaded areas 
are equal to 1/4, 1/2, and 3/4, respectively. Then the respective quartile domains will be 
represented by the intervals on the threshold axis confined between the left and the right 
edges of these areas, and the respective quartile volumes will be the sums of the lengths of 

20 these intervals. 

In Figs. 12 a and 12 b, the quantile densities, volumes, and domains are displayed as 
time dependent quantities computed in a 1-second rectangular sliding window. Panels la 
and 4a of Fig. 12 a show the median densities, computed for the amplitude and the counting 
densities of the signal from Panel la of Fig. 10. Panels lb and 4b of Fig. 12 b show the 

25 respective median densities for the signal from Panel lb of Fig. 10. Panels 2a and 5a of 
Fig. 12 a, and Panels 2b and 5b of Fig. 12 b, show the median volumes of the amplitude 
and the counting densities of the respective signals. As can be seen from these examples, 
both quantile densities and quantile volumes characterize the total width of the densities, 
that is, the total size of high density regions in the threshold space. Panels 3a and 6a of 

30 Fig. 12 a, and Panels 3b and 6b of Fig. 12 b, display the quartile domains, with the median 
domain shaded by the gray color, the 5 = 3/4 domain shaded by the light gray, and the first 
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quartile domain shaded black. These examples illustrate how the quantile domain reveals 
the location of the high density regions in the threshold space. 

Earlier in this disclosure, we adopted the restricted definition of a "phase space" as the 
threshold space of the values of the variable, complemented by the threshold space of the 

5 first time derivative of this variable. Thus for a scalar variable the modulated threshold 
densities in such phase space are the two-variate densities. The introduction of the phase 
space densities expands the applicability of the density analysis, and allows more detailed 
study of the changes in the variables. For example, Fig. 13 provides an illustration of the 
sensitivity of the phase space threshold densities to the signal's shape. The first column 

10 of the panels in the figure shows the fragments of three different signals in rectangular 
windows. The second column of the panels shows the phase space amplitude densities 

b{D^, Di, t) = {do^J'^D. - x{s)] dv.Ti^D^ - x{s)\ )^ , (60) 
and the third column displays the phase space counting densities 

r{D,,D,, t) = ^ , , ==T ^ . (61) 

15 This figure also illustrates that while the amplitude density is indicative of the "occupancy" 
(that is, the time the variable occupies the infinitesimal volume in the phase space), the 
counting density reveals the "traffic" in the phase space, that is, it is indicative of the rates 
of visiting a small volume in the phase space. 

The example in Fig. 14 shows time dependent phase space amplitude densities computed 

20 according to Eq. (60) in a 1-second rectangular moving time window for two computer 
generated non-stationary signals shown in Panels la and lb of Fig. 10. The figure plots the 
level lines of the phase space amplitude densities (Panels la and 2a), at times indicated by 
the time ticks. Panels lb and 2b show the time slices of these densities at time t = to. 
Fig. 15 shows time dependent phase space counting rates 



(62) 



computed in a 1-second rectangular moving window for the two signals shown in Panels la 
and lb of Fig. 10. The figure plots the level lines of the phase space counting rates (Panels la 
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and 2a) at times indicated by the time ticks. Panels lb and 2b show the time slices of these 
rates at time t = to. As was discussed in the explanation of Fig. 10, the general shape of the 
pulses around this time for both signals is similar. Thus the difference in the phase space 
crossing rates apparent from these time slices results mostly from the small differences in 

5 the shape of the "tops" of these pulses. 

Figs. 16 and 17 display the boundaries of the median domains for the phase space 
amplitude and counting densities, respectively. The upper panels of these figures show the 
respective boundaries for the signal of Panel la of Fig. 10, and the lower panels show the 
median domain boundaries for the signal of Panel lb of Fig. 10. 

10 The examples in Figs. 9 through 17 illustrate the usefulness of AVATAR for visual 
assessment of various features of a signal, and of the evolution of these features in time. 
Although Figs. 13 through 17 deal with the phase space densities of a scalar variable, it 
should be obvious that densities of any two-dimensional variable can be treated in the same 
way. For instance, the same technique will apply for describing the time evolution of the 

15 population in a geographic region (amplitude density) and for mapping out the routes of 
migration of this population (counting density). Fig. 16, for example, can represent the 
time evolution of a quantile domain of the population of a biologic species, that is, the 
locations of the largest number of specimens per unit area in a region. Then Fig. 17 will 
represent the time evolution of the respective quantile domain of the traffic of the species, 

20 that is, the regions of the most active movement of the species. 

11.1 ELIMINATING THE DIGITIZATION-COMPUTATION STEPS IN 
DISPLAYING DENSITIES AND THEIR TIME EVOLUTION: DIRECT 
MEASUREMENT BY CONTINUOUS ACTION MACHINES 

Since many physical sensors have input-output characteristics equivalent to those of the 
probes in this disclosure, there is a large number of feasible physical embodiments of 
AVATAR for displaying the modulated threshold densities and their time evolution by 
25 continuous action machines. In this subsection, we provide several illustrative examples 
of such embodiments. The underlying motivation behind constructing an analog machine 
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directly displaying the densities is schematically stated in Fig. 18. 

As the first example, Fig. 19 outlines a conceptual schematic of a simple device for 
displaying time dependent amplitude densities of a single scalar signal. An electron gun in 
a cathode-ray tube produces a beam of fast electrons. The tube contains a pair of vertical 

5 deflecting plates. By feeding a voltage to this pair of plates, we can produce a proportional 
displacement of the electron beam in the vertical direction. The screen of the tube is coated 
with luminophor with the afterglow half-time T1/2 = Tln2. We assume that the brightness 
of the luminescent spot on the screen is proportional to the intensity of the electron beam, 
and is described by dYJ^AV^X) when the voltage across the deflecting plates is zero. Then 

10 the brightness of the displayed band on the screen, at any given time, will correspond to 
the amplitude density of the input signal x{t), computed in the exponential moving window 
of time {RCiQ ) with the time constant T. This band can then be projected on a screen 
by, for example, a (concave) mirror M. By rotating this mirror, we can display the time 
evolution of the amplitude density of x{t). If we now modulate the intensity of the electron 

15 beam by the signal K{t), then the brightness of the displayed picture will be proportional 
to {K dyJ^AYiy - x))t- For example, when K{t) = \x{t)\, the screen will display the 
threshold crossing rates. A simple conceptual schematic of such a device for displaying 
time dependent threshold crossing rates of a signal is illustrated in Fig. 20. Note that by 
displaying only the lines of equal intensity, or by thresholding the intensity, we will reveal 

20 the boundaries of the respective quantile domains. 

Fig. 21 provides an illustration for possible hardware implementation of a device for 
displaying time slices of the phase space amplitude densities. An electron gun in a cathode- 
ray tube of an oscilloscope produces a beam of fast electrons. The tube contains two pairs of 
mutually perpendicular deflecting plates. By feeding a voltage to any pair of plates, we can 

25 produce a proportional displacement of the electron beam in a direction normal to the given 
plates. The screen of the tube is coated with luminophor with the afterglow half-time Tif2 = 
Tin 2. We assume that the brightness of the luminescent spot on the screen is proportional 
to the intensity of the electron beam, and is described by dx^^Axi^) dyJ^AyiX) when 
the voltage across the deflecting plates is zero. If the input signals are x{t) and x{t), 

30 respectively, then the displayed picture on the screen, at any given time, will correspond 
to the phase space amplitude density of the input signal x{t), computed in the exponential 
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moving window of time {RCw ) with the time constant T. Thus, the screen will display 
figures similar to those shown in the second column of the panels in Fig. 13, and in Panels 
lb and 2b of Fig. 14. If we now modulate the intensity of the electron beam by the signal 
K{t), then the brightness of the displayed picture will be proportional to {K dx^Ax{X — 

5 x) dyT^viY - x))t. For example, when K{t) = y/{x /\Yf + [xAXf, the screen will 
display the time slices of the phase space threshold crossing rates of the input signal x[t)^ 
computed in the exponential moving window of time (RCio ) with the time constant T. 
Thus, the screen will display figures similar to those shown in the third column of the 
panels in Fig. 13 and in Panels lb and 2b of Fig. 15. A simple conceptual schematic of such 

10 a device for displaying time shoes of the phase space threshold crossing rates is illustrated 
in Fig. 22. 

One skilled in the art will now recognize that the task of displaying the densities can also 
be achieved by a variety of other physical devices. In addition to displaying the modulated 
densities and their time evolution, these devices can also be modified to display the quantile 
15 domain, density, and volume, the means at reference thresholds, various other quantities 
of the MTDs such as their level lines, and to accomplish other tasks of AVATAR. Some of 
the additional embodiments of such devices will be described later in this disclosure. 



12 USING MEANS AT REFERENCE THRESHOLDS FOR 
DETECTION AND QUANTIFICATION OF CHANGES IN 
VARIABLES 

As has been indicated earlier in this disclosure, a comparison of the mean at reference 
20 threshold with the simple time average will indicate the interdependence of the input and 
the reference variables. For the purpose of this disclosure, we will adopt one of m.any 
possible ways to measure such dependence within the quantile domain as follows: 

S,(D,O^ I{™g-;)-^^>^l .[(/.(D-x)),-/,W] . (63) 

where fg{t) is the quantile density defined by Eq. (24), and we will assume that the norm is 
25 computed simply as the distance in the Euclidean space. Eq. (63) represents an estimator 
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of differences in the quantile domain between the mean at reference threshold and the time 
average. 

Fig. 23 displays the values of the estimator Eq{'D,t) of Eq. (63) in q — 9/10 quantile 
domain, computed for the two computer generated nonstationary scalar signals (Panels la 
5 and lb), used in a number of our previous examples. Panels 2a and 2b display the values 
of the estimator for K = and Panels 3a and 3b display these values for K = 



As follows from Eq. (9), the time dependent Modulated Cumulative Threshold Distribution 
(MCTD) C/^(D,t), or Threshold-Time Averaged Cumulative Distribution, is defined as 



where CKCD^t) is the modulated threshold density given by Eq. (52). It is easy to see 
from Eq. (64) and from the definition of the Heaviside unit step function, Eq. (4), that all 
equations for the counting, amplitude, and acceleration densities (for example, Eqs. (54), 
(55), (60), and (61)) are valid for the cumulative distributions as well, provided that the 
15 symbols '6', V, 'c', 'S', and 'S^J^' in those equations are replaced by 'B\ 'i?', 'C, '9\ and 
respectively. 

Note that the transition from the densities to the cumulative distribution functions is 
equivalent to the threshold integration of the former, and thus the principal examples of the 
embodiments for the densities can be easily modified for handling the respective cumulative 
20 distributions. For instance, the embodiments of Figs. 19 and 20 can be easily converted 
to display the cumulative distributions instead of the densities. Then, for example, the 
lines of equal intensity on the screen will correspond to the level lines of these cumulative 
distributions, and thus will be equivalent to the outputs of the respective rank filters. 
Note also that the ability to compute or measure the time dependent threshold densities 
25 and cumulative distributions for a signal gives access to a full range of time dependent 
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10 




{K(S) ntl^AD, [Pi - X,{s)])j, 



(64) 
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statistical measures and estimates, such as different threshold moments (e.g., mean and 
median, skew, kurtosis, and so on) with respect to these distributions. Although the defining 
equations for the threshold densities and cumulative distributions are given for a continuous 
variable or signal, in numerical computations these quantities can be calculated in finite 
5 diff'erences. Clearly, the introduction of threshold averaging alleviates the computational 
problems caused by the singularity of the delta function. 

14 UNIMODAL APPROXIMATIONS FOR IDEAL MODULATED 
DENSITIES AND CUMULATIVE DISTRIBUTIONS 

For time independent thresholds, numerical (or hardware) computation of densities and 
cumulative distributions according to Eqs. (52) and (64) should not cause any difficulties. 

10 However, in the equations for rank normalization, Eqs. (12) and (13), and for filtering, 
Eqs. (17) and (18), densities and cumulative distributions appear with thresholds dependent 
on time, and their evaluation may present a significant computational challenge. 

This challenge is greatly reduced when the thresholds vary slowly, which is usually 
the case in rank filtering. In such cases, it might be sufficient to approximate (replace) 

15 (/r [D(t) - x(s)])t by the first term in its Taylor expansion, namely by (/r [D(s) - x(s)])r, 
and higher order terms can be retained when necessary. This approximation will be our 
prime choice in most embodiments of AVATAR discussed further in this disclosure. For rank 
normalization, however, this approach is only adequate when certain relations between the 
input and reference signals are held, and thus different approximations should be developed. 

20 Let us first develop a unimodal approximation for the ideal modulated density function 
({•^)t) ^ - ^))t^ that is, the density function resulting from the measurements 

by an ideal probe. Although we will later present more accurate approximations for rank 
normalization, the unimodal approximation has certain merits of its own, e.g., simphcity 
of implementation and analytical treatment of the results. To develop such unimodal ap- 

25 proximation, we can use, for example, the integral representation of the delta function as 
(Arfken, 1985, p.799, for example) 




(65) 



67 



Atty. Docket No. : 1B252.2 
First Named Inventor: Nikitin 
US Express. EL388608922US 



The expression for the density can thus be rewritten as 

^<ifMI>-x)>, = i£d.e-if^, (66) 
and the time averages in the right-hand side of Eq. (66) can be expressed as 

(^^j^ = e " (k)t ^ = e^(") . (67) 

5 It is worth mentioning that the cumulant function A{u) = ln(^{K):^^ (^K e'''*^'^^^ corre- 
sponds to the thermodynamic characteristic function (free energy divided by kT) in statis- 
tical mechanics (Kubo et al, 1995, for example). We now notice that the real part of A(n) 
has a global maximum at n = 0. Therefore, the main contribution of the integrand to the 
integral in Eq. (66) may come from the region around u = 0 (see, for example, Erdelyi, 

10 1956, Copson, 1967, or Arfken, 1985). Thus we can expand A{u) in Taylor series around 
w = 0 as 

where 

Kn = {Kx'')j,. (69) 
15 Truncating this expansion after the quadratic term, we substitute the result in Eq. (66) 
and easily arrive at the following expression: 

1 exp[-|f^l 

where Kno = Kn Kq^. 

Unimodal approximations for higher-dimensional densities can be developed in a similar 
20 manner. For example, 

{6 {D^ - x) S {Dy - y))^ = ^ du e™^^ dv e^^^ e^^^"'^) , (71) 

where 



A{u,v) = ln(e-'''^-'^^)^ « ~i[^/(2;)-f - 



;^2[(a;2) - {xY] -f 2uv[{xy) - {x){y)] + v^[{y^) - [yf]} . (72) 
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This results in 




(73) 



The approximations of Eqs. (70) and (73) are somewhat accurate, for example, on a 
large time scale, when the signals x(t) and y{t) represent responses of (linear) detector 
systems to trains of pulses with high incoming rates (Nikitin et al, 1998; Nikitin, 1998, for 
example), Poisson distributed in time. Fig. 24 illustrates the adequacy of the approximation 
of Eq. (73) for two-dimensional amplitude densities of such signals. The top panel of this 
figure shows the signals x{t) and y{t) = x{t) in a rectangular window of duration T, where 
x{t) is the response of an RCu filter to a pulse train with an energy spectrum consisting 
of two lines of equal intensity with one having twice the energy of the other. The lower 
left panel shows the measured density (Od^J^ad^ {D^ - x) do^J't^Dy {Py - y))^ (with small 
/\Dr^ and ^Dy), and the lower right panel shows the density computed through Eq. (73). 

More generally, Eqs. (70) and (73) are accurate when x{t) and y{t) (on a time scale T) 
are sums of a large number of mutually uncorrected signals. Obviously, when x{t) and y{t) 
are (nonlinear) functions of such sums, simple transformations of variables can be applied 
to modify Eqs. (70) and (73). For example, if Eq. (70) is adequate for the signal x{t) and, 
in addition, KIq < K20, then the following approximations are also adequate: 



j^{K5{D-z))^^ 



0{D) 



(74) 



V27r D K, 



for the signal z{t) 



x^{t), and 




(75) 



for z{t) = \x{t)\. 

When the signal x{t) is a sum of two uncorrelated signals Xi{t) and X2{t), 




(76) 



69 



us Express. EL388608922US 

and thus the resulting density is just the convolution of the densities of the components 
xi{t) and X2(t). This approach can be used, for example, when the signal is a sum of a 
random noise component and a deterministic component with known density functions. 
Fig. 25 illustrates this for the noisy signal xi = sm(t). The signals Xi{t), X2{t), and 
xi{t) + X2{t) are shown in the left column of the panels, and the respective panels in the 
right column show the amplitude densities. The signal X2{t) is random (non-Gaussian) 
noise. The amplitude density of the sinusoid Xiit) in a square window of length T = 27rw 
is computed as 



{6[D-x, {s)] )t = ^ fj^^^t S [D - cos(i)] . 
2^ S Jo 



|sin(^;i.)| TTsin [arccos(L>)] ' 



tk = (-1)^ arccos(i?) + | [2A; + 1 - (-l)'^] . (78) 

In the lower right panel, the measured density of the combined signal is shown by the solid 
line, and the density computed as the convolution of the densities bi{D) and b2{D) is shown 
by the dashed line. 

Substitution of Eq. (70) into Eq. (64) leads to the approximation for the time dependent 
cumulative distribution as follows: 



(79) 



where erfc(a:) is the complementary error function (Abramowitz and Stegun, 1964, for 
example) . 

The unimodal approximations of this section are of hmited usage by themselves, since 
they are adequate only for special types of a signal on a relatively large time scale. For 
example, the approximation of Eq. (70) is generally a poor approximation, since every 
extremum of the signal x{t) inside of the moving window may produce a singularity in 
the density function ((i^)^)"^ {K S {D - x))^,. This can be clearly seen from Eq. (40), 
which shows that a modulated density might be singular whenever \x(ti)\ = 0, with the 
exception of the counting density. For the latter, the approximation of Eq. (70) might be an 
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adequate choice for some usages. This is illustrated by Fig. 26, which shows the amplitude 
{6 (D - x))^ and the counting {{\x\)^)-^ {\x\ 5 [D ~ x))^ densities of the fragment of the 
signal displayed in the upper panel. One can see that the Gaussian unimodal approximation 
(dashed lines) is more suitable for the counting density than for the amplitude density. The 
latter is singular at every threshold passing through an extremum of x{t). 

However, the unimodal approximations of this section might be a good practical choice 
of approximations for rank normalization, since the main purpose of the latter is just 
providing a "container" in the threshold space, where the reference variable is likely to be 
found. This will be discussed in more detail further in the disclosure. 

15 RANK NORMALIZATION 



Although Eqs. (12) and (13) introduce rank normalization for vector fields, in various 
practical problems this normalization may be more meaningfully applied to scalar variables 
and fields, e.g., to the selected components of vector variables. For the time being, we 
will also consider only single scalar or vector variables rather than the fields. As will be 
illustrated further in the disclosure, the transition to the fields can always be accomplished 
by spatial averaging. 

Let us write a reference threshold distribution, modulated by some associated signal K, 

as 

CkA^, t) j^^^ ( K{s) [D - v{s)] )^ , (80) 

where r(t) is the reference signal, i.e., the signal for which the distribution is computed. 
Then 

y{t) = Ck,. m,t] = { K{s) [x(^) - r(5)] )^ (81) 

is a signal, rank normalized with respect to the reference distribution CK,r{^,t). Eq. (81) 
is thus a defining equation for an Analog Rank Normalizer (ARN). In other words, an 
ARN outputs the rank of the value of x(i) with respect to the sample of the values of the 
reference variable v{t). For example, if the reference distribution for normalization of a 
scalar variable x{t) is provided by a Gaussian process with the mean x and the variance 
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, then 



y{t) = - erfc 



(82) 



where we have neglected the threshold averaging, that is, we assumed that J^ad(x) — 
Using Gaussian approximation for modulated threshold densities of a scalar reference r(i), 
Eq. (70), we can write the following expression for the Gaussian normalization: 



y{t) 



erfc 



(83) 



(84) 



10 Let us now observe that 

yit) = |_"d"D 5 [D - x(0] CkA^, t) ^ |_~d"D /r [D - x(i)] CK,r(D, t) , (85) 

where we replaced the Dirac (5-function c^(x) by the response of a probe /r(x) with a small 
width parameter R. Thus, for a small width parameter R, the practical embodiment of 
ARN reads as follows: 

y{t) = ( /r [D - x(i)] CkA^, t) Z , (86) 
where we used the shorthand notation for the threshold integral as 

= (87) 

Fig. 30 shows a simplified flowchart of an analog rank normalizer according to Eq. (86). 
For the purpose of "confinement" of the variable in the threshold space of the reference 
20 variable, rank normalization can be defined by means of a discriminator with arbitrary 
input-output response, namely by 

y{t) = ^R.(o [D,(t) - x(t)] = n J^AD^^^^t) [DrAt) - x,{t)] , (88) 

where Dr(t) is indicative of the central tendency (such as mean or median) of the refer- 
ence density CK,r(D,t), and Rr(t) is indicative of the width (such as standard or absolute 
25 deviation, or such as FWHM) of the reference density CK,r{'D,t). For instance, using the 
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mean and standard deviation of the reference distribution as the displacement and width 
parameters, respectively, for a scalar variable we have 

^(^) = ^y' 2iK,o-Kf,) [^10 W - X{t)] . (89) 

Eq. (89) thus performs rank normalization as transformation of the input signal by a dis- 
criminator with the width parameter ^2[K2oit) - K^o{t)], and the displacement parameter 
Kw{t). 

Rank normalization can also be accomplished through evaluating the integral of Eq. (81) 
by straightforward means. For instance, for a causal time weighting function, 

(^R [x(t) - r(5)])y = I'ds hit - s) [x(i) - r{s)] = 

= ^°ds h{s) ^R [x{t) - r{t - s)] , (90) 
which for a rectangular time function of duration T leads to 

(J^R[x(t)-r(5)])y = ^l^dsJ^^[^{t)-r{t-s)]^ 

~ TfY:^^[^it)-At-nAt)] , (91) 

where At = T/N. 

One of the main usages of the rank normalization is as part of preprocessing of the 
input variable, where under preprocessing we understand a series of steps (e.g., smoothing) 
in the analysis prior to applying other transformations such as MTD. Since in AVATAR the 
extent of the threshold space is determined by the reference variable, the rank normalization 
allows us to adjust the resolution of the acquisition system according to the changes in the 
threshold space, as the reference variable changes in time. Such adjustment of the resolution 
is the key to a high precision of analog processing. In the next section, we provide several 
examples of the usage of rank normalization. 
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16 USING RANK NORMALIZATION FOR COMPARISON OF 
VARIABLES AND FOR DETECTION AND QUANTIFICATION 
OF CHANGES IN VARIABLES 

The output of a rank normalizer represents the rank of the test signal with respect to the 
reference distribution. Thus comparison of the outputs of differently normalized test signals 
constitutes comparison of the reference distributions. Various practical tasks will dictate 
5 different implementations of such comparison. Let us consider several simple examples of 
using rank normalization for comparison of reference distributions. 

First, let us define a simple estimator Qabit] q) of differences between the distributions 
Ca{D,t) and Cb{D,t) as 

Qab{t;q)=^a[y,{t),t] 
CaMt),t]^q 

10 In Eq. (92), yq(t) is the level line {qth quantile) of the distribution Ca{D,t). Clearly, when 
Ca{D, t) and Cb{D, t) are identical, the value of Qab{t] q) equals the quantile value q. Fig. 27 
provides an example of the usage of the estimator given by Eq. (92) for quantification of 
changes in a signal. In this example, the signals are shown in Panels la and lb. The distri- 
butions Ca{D,t) are computed in a 1-second rectangular moving window as the amplitude 

15 (for Panels 2a and 2b) and counting (for Panels 3a and 3b) cumulative distributions. Thus 
yq{t) are the outputs of the respective rank filters for these distributions. The estimators 
Qab{t; q) are computed as the outputs of the Gaussian normalizer of Eq. (83). The values of 
these outputs for different quartile values are plotted by the gray (for q = 1/2), black (for 
q = 1/4), and light gray (for q = 3/4). In this example, the estimator Qab(t;q) quantifies 

20 the deviations of Ca{D, t) from the respective normal distributions. 

Since the shape of the signal has different effects on the amplitude, counting, and 
acceleration distributions, comparison of the signal normalized with respect to these 
distributions can allow us to distinguish between different pulse shapes. Fig. 28 provides 
a simplified example of the usage of rank normalization for such discrimination between 
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different pulse shapes of a variable. Panel I shows the input signal consisting of three 
different stretches, 1 through 3, corresponding to the variables shown in Fig. 9. Panel II 
displays the difference between C|^|_^(a;,t) and C'l^^{x,t), where the superscripts ho denote 
the particular choice of the time weighting function as an Rdo filter, and the reference 
5 signal r is a Gaussian process with the mean Kiq and the variance 7^20 - K^q, where 
Knm are computed for the input signal x{t). This difference is zero for the first stretch 
of the input signal, since for this stretch the amplitude and the counting densities are 
identical (see Fig. 9). Panel III displays the difference between C|^|^(a7,t) and C^°,{x,t). 
This difference is zero for the second stretch of the input signal, since for this stretch the 
10 amplitude and the acceleration densities are identical (see Fig. 9). The distance between 
the time ticks is equal to the constant T of the time filter. Fig. 31 shows a simplified 
flowchart of a device for comparison of two signals. In order to reproduce the results shown 
in Fig. 28, Panels II and III, the specifications of the device are as follows: 

xi{t) = X2{t) = n{t) = r2{t) = x{t)] 

Ki{i) = constant ; 

< K2{t) = \xit)\ for Panel II, K2{t) = \xit)\ for Panel III; (93) 

^Ao(a;) = ierfc(-^) ; 
_ h{t) = ho(t) {RCio filter) . 

15 

Fig. 29 provides an additional example of sensitivity of the difference between two rank 
normalized signals to the nature of the reference distributions. Panel I shows the input 
signal, and Panel II displays the amplitude density computed in an RCw window with the 
time constant T, equal to the distance between the time ticks. Panel III plots the difference 

20 between CfP^^^ {x, t) and C^^^ {x,t), where the reference signal r is a Gaussian process with the 
mean Kiq and the variance K20 - K^q, computed for the input signal x{t). The magnitude 
of this difference increases due to broadening of the amplitude density while the counting 
density remains unchanged. 

Even though various practical tasks will dictate different implementations of comparison 

25 of variables through rank normalization, the simple examples provided above illustrate 
that sensitivity of the difference between two rank normalized signals to the nature of the 
reference distributions provides a useful tool for development of such implementations. 
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16.1 ESTIMATORS OF DIFFERENCES BETWEEN TWO CUMULATIVE 
DISTRIBUTIONS 



Since AVATAR transforms a variable into density and cumulative distribution functions, an 
expert in the art would recognize that any of the standard techniques for statistical analysis 
of data and for comparison of distributions/densities can also be applied to time dependent 
5 quantification of the signal and its changes. For example, one can use an estimator of the 
differences between the two distributions Ca and Cb as 

nt) = ^m)] = ^[Km , (94) 

where <^ is some function, and A„6 is the statistic of a type 

AatW =/_^d"x/i(x)i^[a(x,t),a(x,t),c,(x,t),c6(x,t)] , (95) 

10 where /i is a truncation function, and H is some score function. For convenience, we 
shall call Ca the test distribution, and Cb the reference distribution. The statistics of 
Eq. (95) thus quantify the differences and changes in the input signal with respect to 
the reference signal. One can think of any number of statistics to measure the overall 
difference between two cumulative distribution functions (Press et al., 1992, for example). 

15 For instance, one would be the absolute value of the volume between them. Another 
could be their integrated mean square difference. There are many standard measures of 
such difference, e.g., Kolmogorov-Smirnov (Kac et al., 1955), Cramer-von Mtses (Darling, 
1957), Anderson-Darling, or Kuiper's statistics, to name just a few (Press et al., 1992, for 
example). For example, one can use the statistics of a Cramer-von Mises type (Darhng, 

20 1957, for example) 

Kt, = I^JCaix) wlCa(x)] W[Ca{x) - C^ix)] , (96) 

where w is a truncating function, and W is some score function. Generalization of Eq. (96) 
to many dimensions is straightforward as 



76 



us Express EL38S608922US 




(97) 



where Ca(x, t) is the density function. 

If Ca is computed for the same signal as C5, then the difference between Ca and Cb is 
due to either the different nature of Ca and (e.g., one is the amplitude, and the other 
5 is counting distribution), or to the difference in the discriminators and time windows used 
for computing Ca and Cb- 

Notice that rank normalization is just a particular special case of the estimator of 
Eq. (97), when w{x) = 1, W{x) = 1/2 — x, and the test density is the instantaneous 
density, Ca(D, i) = 5 [D - x(i)], and thus 



One should also notice that an estimator of differences between two distributions of the type 
of Eq. (94) can be computed as a time average when the distribution functions are replaced 
by the respective rank normalized signals. A simplified example of such an estimator is 
shown in Fig. 32. This figure plots the time averages of the absolute values of the differences, 
15 {\Ck°x{^, t) - C'^U^, ^)I)t, for K = \x\ and K = for the signal shown above the panel 
of the figure. The distance between the time ticks is equal to the time constants T of the 
time filtering windows. 



Selectivity of comparison of the amplitude and counting densities of a scalar variable can 
20 be greatly increased by comparing these densities in the phase space of this variable. Here 
under phase space we understand simply the two-dimensional threshold space for the values 
of the variable as well as the values of its first time derivative. Fig. 33 illustrates sensitivity 
of the amplitude and counting phase space densities to differences in the signal's wave form. 
The panels in the left column show the sound signals for several letters of the alphabet. The 
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16.2 SPEECH RECOGNITION 
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top signals in the individual panels are the original input signals. The normalized input 
signals and their normalized first derivatives, respectively, are plotted below the original 
input signals. The middle column of the panels shows the amplitude, and the right column 
the counting densities of these pairs of normalized signals. Notice that rank normalization 
5 of the components of the signal allows us to more efficiently utilize the threshold space, 
and thus to increase precision of analog processing. Rank normalization also alleviates the 
dependence of the phase space densities on the magnitude of these components. 

Fig. 34 illustrates how statistics of Eq. (97) can be used in combination with rank 
normalization to address the speech recognition problem. Panel I: The original speech 
10 signal "Phase Space" is shown in the top of the panel. This signal is normalized with 
respect to a Gaussian process with the mean and variance of the original signal in a moving 
rectangular window of 45 ms, and the result is plotted just below the original signal. The 
bottom of the panel shows the time derivative of the speech signal, normalized the same 
way. Panel II: Time slices of the threshold density c{D^,Dyj), where x and y are the 
15 normalized original signal and its normalized derivative, respectively, and c{D:,,Dy,t) is 
their amplitude density in the time window 45 ms. The slices are taken approximately 
through the middles of the phonemes. Panel III: Time slices of the cumulative distribution 
C{D^, Dy, t), where x and y are the normalized original signal and its normalized derivative, 
respectively, and C{D^,Dy,t) is their distribution in the time window 45 ms. The slices 
20 are taken approximately through the middles of the phonemes. Panel IV: The value of the 
estimator of a type of Eq. (94), where the reference distribution is taken as the average 
distribution computed in the neighborhood of the phonemes "a" . The larger values of the 
estimator indicate a greater similarity between the signals. 

Employing different variables for analysis, different time weighting windows (i.e., of 
25 different shape and duration), different types of reference distributions for normahzation 
(i.e., Gaussian or of different random or deterministic signal), different functions for spatial 
averaging, different type estimators of the differences between the two distributions (i.e., 
different truncating functions h{x), and different functions H{x) in Eq. (95)), and so on, 
one can reach any desired compromise between robustness and selectivity in identifying 
30 various elements (e.g., phonemes or syllables) in a speech signal. 

As has been discussed previously, the embodiments of AVATAR allow their implemen- 
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tation by continuous action machines. For example, Fig. 35 outlines an approach one may 
take for eliminating the digitization-computation steps in the analysis and for direct imple- 
mentation of speech recognition in an optical device. Imagine that we can modulate the 
intensity of a beam of Hght in its cross section by J^ax(^ - x) J^c\y{Y - y), where x{t) 

5 and y{t) are proportional to the components of the input signal. This can be done, for 
example, by moving an optical filter with the absorption characteristic ^ax{^) ^Ar{^) in 
the plane perpendicular to the beam, as illustrated in Panel I of Fig. 35. In the example 
of Fig. 35, the two components of the input signal are taken as squared rank normalized 
original speech signal, and its squared rank normalized first time derivative. A second beam 

10 of Hght, identical to the first one, is modulated by 1 - J^^^xi^ - r) T^yiy - q) (Panel II), 
where r{t) and q{t) are proportional to the components of the reference signal. In the 
example of Fig. 35, the reference signal is taken as the input signal in the neighborhood 
of the phonemes "a". These two light beams are projected through an (optional) optical 
filter with the absorption characteristic h{X,Y) (Panel III) onto the window of a photo- 

15 multiplier, coated with luminophor with the afterglow time T (Panel IV). We assume that 
the photomultiplier is sensitive only to the hght emitted by the luminophor. Therefore, the 
anode current of the photomultiplier will be proportional to 

A{t)^ JJdXdYh{X,Y) {l~J^Ax{X-r)TAY{Y-q)+rAx{X -x)y^AY{y-y))T , (99) 
and the variance of the anode current on a time scale T will be proportional to {(A — 1)^)t- 
20 Thus any measurement of the variance of the anode current will be equivalent to compu- 
tation of the estimator of Eq. (94), as illustrated in Panel V of Fig. 35. In this example, 
the larger variance corresponds to a greater similarity between the test and the reference 
signals. 

16.3 PROBABILISTIC COMPARISON OF AMPLITUDES 

25 As an additional example of usefulness of the rank normalization for comparison of signals, 
consider the following probabilistic interpretation of such comparison. For a nonnegative 
time weighting function h{t), J^^dth{t) = 1, we define x{s) to be a value drawn from 
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a sequence x{t), provided that s is a random variable with the density function h{t - s). 
When h{t — s) = gi{t — s), then x{s) is a value drawn from the first sequence, and when 
h{i — s] = g2{t — s), x{s) is a value drawn from the second sequence. Now consider the 
following problem: What is the (time dependent) probability that a value drawn from the 
5 first sequence is q times larger than the one drawn from the second sequence? Clearly, this 
probability can be written as 

= l^dy t) C!% (^^, t^ . (100) 

Substituting the expression for Ci\^{y,t) (without spatial averaging, for simplicity), 

c{\Ay.i) = f^Js g^{t - s) 5[y - x{s)] , (101) 

10 into Eq. (100) leads to 

P,{t) = /_^d. g,{t - s) Cl% = {ci% ^1 Y ■ (102) 

In Eq. (102), Cf_^^ [^'*] is the result of the rank normalization of the sequence x{t)/q, 
with respect to x{t), in the second time window g2{t-s). Pg{t) is thus a simple time average 
(with the first time weighting function) of this output. 

15 The probabilistic estimator Pg{t) can be used in various practical problems involving 
comparison of variables. For example, one of the estimators of changes in a nonstationary 
sequence of data can be the ratio of the medians of the amplitudes of the data sequence, 
evaluated in moving time windows of sufficiently different lengths. This approach is used, for 
instance, in the method by Dorfmeister et al. for detection of epileptic seizures (Dorfmeis- 

20 ter et al., 1999, for example). In this method, the onset of a seizure is detected when the 
ratio Xi{t)lx2{t), where Xi{t) is a median of the squared signal in the ith window, exceeds 
a predetermined threshold q. It is easy to show that the inequality Xi{t)/x2{t) > q corre- 
sponds to the condition 

P<i{t)>l, (103) 
25 where Pg{t) is given by Eq. (102), in which x{t) is the squared input signal. Thus the 
computational cost of the seizure detection algorithm by Dorfmeister et al can be greatly 
reduced, and this algorithm can be easily implemented in an analog device. 
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17 ANALOG RANK FILTERS (ARFs) 

Id some applications, one might be interested in knowing the quantile function for the 
signal, that is, in knowing the value of Dg(t) such that Cic(D^,t) = q = constant, where 
Cif(D,i) is a cumulative distribution function. Thus D^(i) is an output of an analog 
5 rank (also order statistic, or quantile) filter. For example, Di/2(^) is the output of an 
analog median filter. Notice that, since the partial derivatives of CK(D,t) with respect to 
thresholds are nonnegative, CK(Dq,t) = q describes a simple open surface in the threshold 
space. 

For amplitudes of a scalar signal, or ensemble of scalar signals, numerical rank filtering is 

10 a well-known tool in digital signal processing. It is a computationally expensive operation, 
even for the simple case of a rectangular moving window. First, it requires knowing, at any 
given time, the values of latest data points, where A'' is the length of the moving window, 
times the number of signals in the ensemble, and the numerical and chronological order of 
these data points. In addition to the computational difficulties due to the digital nature 

1-5 of their definition, this memory requirement is another serious obstruction of an analog 
implementation of rank filters, especially for time weighting windows of infinite duration. 
Another computational burden on different (numerical) rank filtering algorithms results 
from the necessity to update the numerically ordered list, i.e., to conduct a search. When 
the sampling rate is high, A'' can be a very large number, and numerical rank filtering, as 

20 well as any analog implementation of digital algorithms, becomes impractical. 

In AVATAR, the output of a rank filter for a scalar variable is defined as the thresh- 
old coordinate of a level line of the cumulative distribution function. Since the partial 
derivatives of the latter with respect to both threshold and time are enabled by definition, 
the transition from an implicit CK{Dg,t) = q form to the explicit Dg = D^{t) can be 

25 made, for example, through the diflferential equation given by Eq. (11) (see Bronshtein and 
Semendiaev, 1986, p. 405, Eq. (4.50), for example). The differentiability with respect to 
threshold also enables an explicit expression for the output of a rank filter, for example, 
through Eq. (23). Thus AVATAR enables implementation of order statistic analysis in ana- 
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log devices, and offers significant improvement in computational elSciency of digital order 
statistic processing. 



18 ANALOG RANK FILTERS OF SINGLE SCALAR VARIABLE 



AVATAR enables two principal approaches to analog rank filtering, unavailable in the 
prior art: (1) an explicit analytical expression for the output of an ARF, and (2) a differ- 
ential equation for this output. In this section, we briefly describe these two approaches. 

18.1 Explicit Expression for Output of Analog Rank Filter 

An explicit expression for the output of an analog rank filter can be derived as follows. 
Notice that 



where Dq is a root of the function Ck{D, t) - q, 0 < q < 1. Since, at any given time, there 
is only one such root, we can use Eq. (27) (Rumer and Ryvkin, 1977, p. 543, for example) 
to rewrite Eq. (104) as 



where we replaced the Dirac (5-function S{x) by the response of a probe dqJ^Aq{x), and used 
the shorthand notation of Eq. (87) for the threshold integral. Note that the rank filter 
represented by Eq. (105) can be implemented by continuous action machines as well as by 
numerical computations. 

Fig. 36 illustrates the relationship between the outputs of a rank filter and the level lines 
of the amplitude distribution of a scalar signal. Panel I of the figure shows the input signal 
x{t) on the time-threshold plane. This signal can be viewed as represented by its instan- 
taneous density 5[D-x{t)]. Threshold integration by the discriminator Tad{D) trans- 
forms this instantaneous density into the threshold averaged distribution T^d [D - x{t)] 




(104) 




(105) 
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(Panel II). This distribution is further averaged with respect to time, and the resulting 
distribution B{D,t) = {:F^n [D - x{s)])j, is shown in Panel III. The quartile level lines are 
computed as the outputs of the rank filter given by Eq. (105), and are plotted in the same 
panel. Panel IV shows the input signal x{t), the level lines of the amplitude distribution 

5 for = 1/4, 1/2, and 3/4 (gray lines), and the outputs of a digital rank order filter (black 
lines). It can be seen from this panel that the outputs of the respective analog and digital 
rank filters are within the width parameter AD of the discriminator. 

Fig. 37 repeats the example of Fig. 36 for the respective analog and digital median 
filters for the discrete input signals. The instantaneous density of a discrete signal can be 

10 represented by S [D - x{t)]'£iS{t - t^), as shown in Panel I. Panel II shows the thresh- 
old averaged distribution J^ad [D - x(t)] Y.i5{t- t^), and Panel III of the figure compares 
the level line B{D,t) = (^ac [I? - ^(s)] Ei <5(s - ^«))r = 1/2 (solid black line) with the 
respective output of a digital median filter (white dots). 

Fig. 38 shows a simplified schematic of a device for analog rank filtering according to 

15 Eq. (105). 

18.2 Differential Equation for Output of Analog Rank Filter 



Substituting the expression for the modulated cumulative threshold distribution function, 
Eq. (64), into Eq. (18), we arrive at the differential equation for an analog rank filter of a 
single scalar variable as follows: 

p _ d,CK{D„t) _ W^[g-C>(g„i)] 

c,(D,.t) - (K%dl,(D„t) ' (""^ 

where the dots over Dq and h denote the time derivatives, and we used the fact that 
^KiDg^t) = q. In Eq. (106) we used the superscripts h and h to indicate the particular 
choice of the time weighting functions in the time integrals. 

Notice that if h(t) is a time impulse response of an analog filter, then Eq. (106) can 
be solved in an analog circuit, provided that we have the means of evaluating C^{Dq,t) 
and c\{Dq,t). The example in Fig. 39 shows a simplified schematic of such a device for 
analog rank filtering. Module I of the device outputs the signal {K)i^ [q - C^{Dq,t)], and 
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Module II estimates {K)j'c\{Dq,t). The outputs of Modules I and II are divided to form 
Dq{t), which is integrated to produce the output of the filter Dq{t). 
Notice also that in the absence of time averaging 

CK{Dg,t) = c'j,{Dq,t) = OdJ'ad [Dq[t) - x{t)] , (107) 

5 and 

dtCK{Dq,t) = dtC'K{Dq,t) = Ci{Dq,t) = -x{t) ^p^AD [Dq{t) - x{t)] , (108) 

and thus Eq. (106) is still valid, although it leads to the trivial result Dq{t) = x(t) for any 
value of q. 



20 



19 RCin ANALOG RANK FILTERS 

Since RC\n time impulse response functions commonly appear in various technical solu- 
tions, and are easily implementable in analog machines as well as in software, they are a 
natural practical choice for h{t). In addition, the exponential factor in these time weighting 
functions allows us to utilize the fact that (e^)' = e^, and thus to simplify various practical 
embodiments of AVATAR. Substitution of Eqs. (58) and (59) into Eq. (106) leads to the 
following expressions for the analog rank filters: 

for n > 1, and 

for the exponentially forgetting {RCw ) filter. 

20 ADAPTIVE ANALOG RANK FILTERS (AARFs) 



For practical purposes, the averages (/r [D(t) — x(s)])t and (J^r [D(^) — x(s)])t in 
the expressions for analog rank filters can be replaced by (/r [D(s) - x(s)])t and 
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(^R,[D(s) — x(s)])j' , respectively. However, since we allow the variable x(i) to change 
significantly over time, the size of the characteristic volume R needs to be adjusted in ac- 
cordance with these changes, in order to preserve the validity of these approximations. For 
instance, the adaptation scheme can be chosen as 

AD = AD{t) =e + ra{t) , (111) 
where s is the minimal desired absolute resolution, <7^ is the variance of the input signal, 
cr^{t) = K20 — K^Q, and r -C 1 is a small positive number. As a rule of thumb, r should be 
inversely proportional to the characteristic time T of the time weighting function. Other 
adaptation schemes can be used as needed. The preferred generic adaptation should be 
such that the width parameter AD{t) is indicative of variability of Dq{t). For example, 
such adaptation can be as 

AD ^ ADit) =e + r [{dI{s))^ - {D,{s))l] , (112) 

where e is the minimal desired absolute resolution, and r <C 1 is a small positive number. 
Then the equation for adaptive analog rank filters (AARFs) reads as 

^ ^ Q {K{s)i - {K{s) :F^nis) - x(s)])^ ^^^^^ 

{K{s)dr>:F^nis)[D,{s)-x{s)])l 

and in the special (due to its simplicity) case of the exponentially forgetting (.RCio ) filter, 
as 

n - K[,-^.,is,iD,-.)] ^^^^^ 



T{K{s)dDJ'^nis) [D,{s)-x{s)])° 

Fig. 40 shows a simplified diagram of the implementation of Eq. (113) in an analog de- 
vice, with the adaptation according to Eq. (111). Module I takes the outputs of Modules II 
and III as inputs. The output of Module I is also a feedback input of Module IL Module IV 
outputs AD{t), which is used as one of the inputs of Module II (the width parameter of 
the discriminator and the probe) for adaptation. 

When Eq. (113) is used for the implementations of the i?Ci„ time window ARFs, the 
resulting algorithms do not only allow easy analog implementation, but also offer significant 
advantages in numerical computations. Since numerical RC in -filtering requires only n + 1 
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memory registers, an RCin moving window adaptive ARF for an arbitrary associated signal 
K{t), regardless of the value of the time constant T, requires remembering only 6n + 4 data 
points, and only 3n + 1 without adaptation. An easy-to-evaluate threshold test function 
can always be chosen, such as Cauchy. This extremely low computational cost of AARFs 

5 allows their usage {either in analog or digital implementations) in systems and apparatus 
with limited power supplies, such as on spacecraft or in implanted medical devices. 

Fig. 41 compares the quartile outputs (for q = 0.25, 0.5, and 0.75 quantiles) of the 
Cauchy test function RCu AARF for signal amplitudes with the corresponding conven- 
tional square window digital order statistic filter. The outputs of the AARF are shown 

10 by the thick black solid lines, and the respective outputs of the square window order 
statistic filter are shown by the thin black lines. The time constant of the impulse re- 
sponse of the analog filter is T, and the corresponding width of the rectangular window 
is 2aT, where a is the solution of the equation a — ln(l + a) = ln(2). The incoming sig- 
nal is shown by the gray line, and the distance between the time ticks is equal to 2aT. 



20.1 ALTERNATIVE EMBODIMENT OF AARFs 

15 

In some cases, after conducting the time averaging in the equation for the AARF, Eq. (113), 
this expression loses the explicit dependence on the quantile value q. This will always be 
true, for example, for the case of a rectangular time weighting function and a constant K. 
This difficulty can be easily overcome by observing that 

20 Ck{D„ t) = lun^ { Ck{D,, .))^V > (115) 

which leads to the approximate expression for the partial time derivative of CK{Dq,t) as 

dtCK{D„t) « ^ [CK{D„t)-q] . (116) 

Thus Eq. (113) can be written for this special case as 

D - ^ ^^^"^^^ ~ (^^'^ ^^'^^^ ~ "^^"^J)t 

AT (^(Kis) dnT^nis) mt) - ' 

25 where AT is small. 

Note that a numerical algorithm resulting from rewriting Eq. (106) as 



(117) 
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q-CK{D,,t) 



(118) 



AT {ck{D,,s))[ 



will essentially constitute the Newton-Raphson method of finding the root of the function 
CK{Dq, t)-q^Q (Press et al., 1992, for example). 

Fig. 42 compares the quartile outputs (for g = 0.25, 0.5, and 0.75 quantiles) of the 
Cauchy test function square window AARF for signal amplitudes with the corresponding 
conventional square window digital order statistic filter. The outputs of the AARF are 
shown by the black sohd lines, and the respective outputs of the square window order 
statistic filter are shown by the dashed fines. The widths of the time windows are T in all 
cases. The incoming signal is shown by the gray fine, and the distance between the time 
ticks is equal to T. 

The alternative embodiment of the AARF given by Eq. (117) is especially useful when 
an AARF is replacing a conventional square window digital order statistic filter, and thus 
needs to replicate this filter's performance. Another way of emulating a digital rank order 
filter by a means of analog rank selectors will be discussed later in the disclosure. 



In various practical problems, it is often convenient to express the measured variable as 
an ensemble of variables, that is, as 



where n(/i) d/i is the weight of the yuth component of the ensemble such that S^^df^n{fj.) = 
1. For such an entity, the threshold averaged instantaneous density and cumulative distri- 
bution can be written as 



21 DENSITIES AND CUMULATIVE DISTRIBUTIONS FOR 



ENSEMBLES OF VARIABLES 




(119) 




(120) 



and 




(121) 
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respectively. We will use these equations further to develop practical embodiments of analog 
rank selectors. 

Eqs. (120) and (121) lead to the expressions for the respective modulated density and 
cumulative distributions as 

c.(D; t. nM) ^ l_yn(,) <^''t-' j^'" " -"t^)' >^ , (122) 



and 



The definitions of Eqs. (122) and (123) will be used further to develop the AARFs for 
ensembles of variables. 



22 ANALOG RANK SELECTORS (ARSs) 

10 

Let us find a gth quantile of an equally weighted discrete set of numbers {x^}. Substitution 
of Eqs. (122) and (123) (with n(/i) = ^ E^Ii Hl^-i)) into Eq. (22) leads to the embodiment 
of an analog rank selector as 

da ' [1 - (1 - q)e-] E^ doJ^ADix, - x,) ' ^''^^ 
15 Fig. 43 illustrates finding a rank of a discrete set of numbers according to Eq. (124). Five 
numbers are indicated by the dots on the X-axis of the top panel. The solid line shows 
the density resulting from the threshold averaging with a Gaussian test function, and the 
dashed lines indicate the contributions into this density by the individual numbers. The 
solid line in the middle panel plots the cumulative distribution. The crosses indicate Xq{a) 
20 and Ti\D[Xq{a)] at the successive integer values of the parameter a. The bottom panel 
plots the evolution of the value of Xq{a) in relation to the values of Xi. 

If we allow the variables to depend on time, {a;^^} = {x^{t)}., then a convenient choice 
for the parameter is the time itself, and we can use Eqs. (19) through (21) to develop 
(instantaneous) analog rank selectors. Setting 

25 Qrix^t) = P da (l)T{t - a) h{x\ a, n{fj)) = {h{x\ a, n{n)))^ , (125) 
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we can rewrite Eq. (21) as 



(126) 



For example, choosing ^^(i — a) in Eq. (125) as 

M - a) = ^ e^ d{t -a) = ho{t - a) (127) 

leads to the relation 

lim gT{x,t) = lim j _dae''l^b{x- a, n{ix)) = b{x; t, n(^)) . (128) 
Then the equation for a discrete ensemble analog rank selector reads as follows: 

^ - .t/T 'l~Ein^ J'aD [Xg jt) - Xj it)] ^ 

(129) 



T{^in,dnTAD[x,{s)-x,{s)])%' ' 
where T is assumed to be small. For digitally sampled data, T should be several times the 

10 sampling interval At. Fig. 44 provides a simple example of performance of an analog rank 
selector for an ensemble of variables. In Panel I, the solid line shows the 3rd octile of a set 
of four variables {xi{t) through Xi{t), dashed lines), computed according to Eq. (129). In 
Panel II, the solid line shows the median (g = 1/2 in Eq. (129)) of the ensemble. The thick 
dashed line plots the median digitally computed at each sampling time. The time constant 

15 of the analog rank selector is ten times the sampling interval. 

Obviously, when {xi{t)} = {x{t - iAt)} and rn = 1/N, Eq. (129) emulates an N-point 
square window digital order statistic filter. Emulation of digital order statistic filters with 
arbitrary window is done by replacing the (uniform) weights 1/A'' by n^, J2ni = 1. An 
expert in the art will recognize that any digital rank filter in any finite or infinite time 

20 window can be emulated through this technique. Fig. 45 compares quartile outputs (for 
q = 0.25, 0.5, and 0.75 quantiles) of a square window digital order statistic filter (dashed 
lines) with its emulation by the Cauchy test function ARS (solid black lines) . The incoming 
signal is shown by the gray line, and the distance between the time ticks is equal to the 
width of the time window T. 

25 Fig. 46 shows a simplified schematic of a device (according to Eq. (129)) for analog rank 
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selector for three input variables. 



23 ADAPTIVE ANALOG RANK FILTERS FOR ENSEMBLES OF 
VARIABLES 



The equation for AARFs, Eq. (113), can be easily rewritten for an ensemble of variables. 
In particular, for a discrete ensemble we have 
^ ^ ^ Q {E^ K,{s)t - (E. n, K,(s) T^p^,^ [D.js) - x,{s)])l ^^^^^ 

For a continuous ensemble, the summation Ei is simply replaced by the respective integra- 
tion. 

Fig. 47 provides an example of performance of AARFs for ensembles of variables. This 
figure also illustrates the fact that counting densities do not only reveal different features 

10 of the signal than do the amplitude densities, but also respond to different changes in the 
signal. The figure shows the outputs of median AARFs for an ensemble of three variables. 
The input variables are shown by the gray lines. The thicker black lines in Panels I and II 
show the outputs of the median AARFs for amplitudes, and the thinner black lines in both 
panels show the outputs of the median AARFs for counting densities. All AARFs employ 

15 Cauchy test function and i?Cio time averaging. The distance between the time ticks in 
both panels is equal to the time constant of the time filters. 



24 MODULATED THRESHOLD DENSITIES FOR SCALAR FIELDS 

In treatment of the variables which depend, in addition to time, on the spatial coordinates 
a, the threshold and the time averages need to be complemented by the spatial averaging. 
For example, the modulated threshold density of a variable representing a scalar field, 
X = x(a, t), can be written as 

■ ' ' = WmU ■ * ^ 
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Notice that the time and the spatial averaging obviously commute, that is, 

((■■■>0T = (<--->5->i = <--->"K. (132) 

In the next section, we use Eq. (131) to present the analog rank selectors and analog rank 
filters for scalar fields. 

5 Fig, 48 shows a simplified diagram illustrating the transformation of a scalar field into a 

modulated threshold density according to Eq. (131). The sensor (probe) of the acquisition 
system has the input-output characteristic Bd^c^d of a differential discriminator. The width 
of the probe is determined (and may be controlled) by the width, or resolution, parameter 
(S.D. The displacement parameter of the probe D signifies another variable serving as the 

10 unit, or datum. In Fig. 48, the input variable z(x, i) is a scalar field, or a component of an 
ensemble of scalar fields. For example, a monochrome image can be viewed as a scalar field, 
and a truecolor image can be viewed as a discrete ensemble of scalar fields. The output of the 
probe then can be modulated by the variable ii'(x, i), which can be of a different nature than 
the input variable. For example, jFC(x, t) = constant will lead to the MTD as an amplitude 

15 density, and i^(x, t) = |i(x, i)| will lead to the MTD as a counting density/rate. Both the 
modulating variable K and its product with the output of the probe KdjyTixD can then 
be averaged by a convolution with the time and the space weighting functions hil:; T) and 
/(x;R), respectively, leading to the averages {K doT^olD — 2:(r,t)])^'^ and (i^(r,i))y^. 
The result of a division of the latter average by the former will be the modulated threshold 

20 density cx(D; x, t). Notice that all the steps of this transformation can be implemented by 
continuous action machines. 

25 ANALOG RANK SELECTORS AND ANALOG RANK FILTERS 
FOR SCALAR FIELDS 

For a scalar field (n-dimensional surface) z = z{x,t), where x = {xi, . . . ,Xn) is an n- 
dimensional vector, Eq. (129) can be easily re-written as an RCio analog rank selector/ filter 
25 for a scalar field 
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• ^ q-{T^D[z,M~z{vJ)])i ^^^^^ 

where, as before (see Eq. (47)), (• • ~ J^^d"'r fR{x-r) ■ ■ • denotes the spatial averaging 
with some test function /r(x). 

Fig. 49 shows a simpHfied schematic of a device according to Eq. (133) for an analog 
5 rank filter of a discrete monochrome surface with 3x3 spatial averaging. 

The explicit expression for an ARE, Eq. (105), can be easily re- written for scalar field 
variables as 

D,{a,t) = l2dDDcK{D;a,t)d,J'^,[CK(D;a,t)-g] , (134) 
and the differential equation for an adaptive analog rank filter for a scalar field will read as 

10 

D,M = ^ ^^'^ , (135) 

{k{v, s) OdJ'adm [D,{a, s) - x{r, s)])^^^ 

where A is the width parameter of the spatial averaging filter. 



25.1 IMAGE PROCESSING: ARSs AND ARFs FOR TWO-DIMENSIONAL 
DIGITAL SURFACES 

The simple forward Euler method (Press et al., 1992, Chapter 16, for example) is 
15 quite adequate for integration of Eq. (133). Thus a numerical algorithm for analog rank 
processing of a monochrome image given by the matrix Z = Zij(t) can be written as 

Qk = Qk-i + {q- F)/4 

F = Em,n Wmn J'aD [Qk-1 " {Z,^m,j-n) k~l] ^ ^^^^^ 
fk = g+^{k-l I' 
g = Em,n Qd^aD [Q/c-1 " (Z,-mJ_n)fc-l] . 

where Q is the q'Qa. quantile of Z , and Wmn is some (two-dimensional) smoothing filter, 
YLm,n'^mn — 1- Employing the Cauchy test function, we can rewrite the algorithm of 
20 Eq. (136) as 
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Qk = Qk-i + iq- F)/fk 
F = I + i Em,n Wmn arctan 



(137) 



Fig, 50 provides a simple example of filtering out static impulse noise from a monochrome 
image (a photograph of Werner Heisenberg, 1927) according to the algorithm of Eq. (137). 

5 Panel 1 shows the original image Z . Panel 2 shows the image corrupted by a random 
unipolar impulse noise of high magnitude. About 50% of the image is affected. Panel 3a 
shows the initial condition for the filtered image is a plane of constant magnitude. Panels 3b 
through 3g display the snapshots of the filtered image Q (the first decile of the corrupted 
one, g — 1/10) at steps n. 

10 Fig. 51 illustrates filtering out time- varying impulse noise from a monochrome image (a 
photograph of Jules Henri Poincare) using the algorithm of Eq. (137). Panels la through Ic: 
Three consecutive frames of an image corrupted by a random (bipolar) impulse noise of 
high magnitude. About 40% of the image is affected. Panels 2a through 2c: The image 
filtered through a smoothing filter, (Z = Em,n "^mn ^i-m,j-n • Panels 3a through 3c: The 

15 rank filtered image Q (the median, g = 1/2). The smoothing filter in Eq. (137) is the same 
used in Panels 2a through 2c. 

26 MODULATED THRESHOLD DENSITIES FOR VECTOR FIELD 
AND ENSEMBLE OF VECTOR FIELDS 



The equation for the modulated threshold density of a scalar field, Eq. (131), can be easily 
extended for a vector field as 

(jr(r, .)/R[D-x(r, s)])y^ 
(ir(r,.))J;i 

where A is the width parameter of the spatial averaging filter, and for an ensemble of 
vector fields as 



C7r(D; a, t) = ' ' ' '''^''^ , (138) 
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't,a 



(139) 



Fig. 52 shows a simpHfied diagram illustrating the transformation of a vector field into a 
modulated threshold density according to Eq. (139). The sensor (probe) of the acquisition 
system has the input-output characteristic /r^ of a differential discriminator. The width 
of this characteristic is determined (and may be controlled) by the width, or resolution, 
parameter R^. The threshold parameter of the probe D signifies another variable serving 
as the unit, or datum. In Fig. 52, the input variable x^(a, t) is a component of an ensemble 
of a vector field. For example, a truecolor image can be viewed as a continuous 3D vector 
field (with the 2D position vector a). The output of the probe then can be modulated by the 
variable K^{a.,t), which can be of a different nature than the input variable. For example, 
ir^(a, t) = constant will lead to the MTD as an amplitude density, and K^{sl, t) = |x^(a,i)| 
will lead to the MTD as a counting density/rate. Both the modulating variable and its 
product with the output of the probe K^f^^ can then be averaged by a convolution with 
the time and the space weighting functions h{t-T) and /(a; R), respectively, leading to the 
averages (ir^/R,(D - x^))^'^ and {K^)"^^^. The result of a division of the latter average by 
the former will be the modulated threshold density c;rjD;a,i). Notice that all the steps 
of this transformation can be implemented by continuous action machines. 

27 MEAN AT REFERENCE THRESHOLD FOR VECTOR FIELD 

The equation for the mean at reference threshold, Eq. (53), can be easily extended for a 
vector field input variable K(a, s) as 




(140) 
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28 ANALOG FILTERS FOR QUANTILE DENSITY, DOMAIN, AND 
VOLUME (AQDEFs, AQDOFs, and AQVFs) 



Notice that the quaatile density, domain, and volume (Eqs. (24) and (25)) are defined for 
multivariate densities, and thus they are equally applicable to the description of the scalar 
variables and fields as well as to the ensembles of vector fields. The quantile density /^(a, t) 
defined by Eq. (24), and the quantile volume Rg{a.,t) defined by Eq. (25), are both good 
indicators of an overall (threshold) width of the density /;<'(x;a, t). The analog filters for 
these quantities can be developed as follows. 

Let us denote, for notational convenience, the density function /^;:(x;a, i) in Eq. (24) 
as z{yi,t). Notice that {z{r,t))l^ = d"r /i<'(r; a, ^) = 1, and thus Eq. (24) can be 
re-written in terms of a modulated cumulative threshold distribution of a scalar field, 
namely as 

CAzM-,t] = hm MM)^..k(t)-.(r,t)])L^ 

= hm — -j-^ '~ = l~q. (141) 

Keeping T and AD in Eq. (141) small, but finite, allows us to write the equation for an 
Analog Quantile Density Filter (AQDEF) as 

. dtC^[z,{t)^^] _ {iZc 'i^'^^^^s) {1-q- TisD [^g(^) - ^(r, g)]})^ _ 



^AzM^A (r^d"r^(r,s)9z,^ADk(^) -^(r,s)])J, 

(r^d-r^(r,s) {1 - d - Tc.u{z^{s) - z{x,s)\})\ 

~ J; ^ . (142) 

(/-^d«r 2(r, s) dnT^D \zq{s) - z(r, s)])^ 

For the exponentially forgetting (.RCio) time filter, h = ho, Eq. (142) translates into 

i^(i) = rood'r ^(r, s) {1 - g - J^^n [z,{s) - z{r, s)]} ^^^^^ 
T (/f°^d"r^(r,.s)aD^AD[^,(s) - ^(r,s)])^° 

When the difference between the density function z{'D,t) and the quantile density Zg{t) is 
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the argument of a discriminator J^ap, the resulting quantity SgCDid-jt). 



(144) 



can be called a quantile domain factor since the surface SgCD; a, = | confines the regions 
of the threshold space where z{'D,t) > Zg{t). Thus Eq. (144) can be used as the definition 
of an Analog Quantile Domain Filter (AQDOF). Integrating over all threshold space, we 
arrive at the expression for an Analog Quantile Volume Filter (AQVF) as follows: 



Fig. 53 shows a diagram of a process for the transformation of the incoming vector field 
10 x(a, t) into a modulated threshold density Ci^(D; a, t), and the subsequent evaluation of the 
quantile density Zq{t), quantile domain factor 5*5(0; a, t), and the quantile volume Rg{a.,t) 
of this density. 

Figs. 54 a and 54 b show the median densities and volumes computed for the amplitude 
and counting densities of the two signals used in several previous examples (see, for example, 
15 Figs. (10) through (12b), (14) through (17), (23), and (27)). These figures compare the 
median densities and volumes computed directly from the definitions (Eqs. (24) and (25), 
gray fines) with those computed through Eqs. (143) and (145) (black lines). Panels la, 2a, 
lb, and 2b relate to the amplitude densities, and Panels 3a, 4a, 3b, and 4b relate to the 
counting densities. 




(145) 
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29 SOME ADDITIONAL EXAMPLES OF PERFORMANCE OF ARSs 
AND AARFs 

29.1 Comparison of Outputs of RCu AARFs for Accelerations, 
Threshold Crossings, and Amplitudes 

Fig. 55 shows the quartile outputs (for q — 1/A through 3/4 quantiles) of the RCw Cauchy 
test function AARFs for the signal amphtudes (Panel I), threshold crossing rates (Panel II), 

5 and threshold crossing accelerations (Panel III). The signal consists of three different 
stretches, 1 through 3, corresponding to the signals shown in Fig. 9. In Panels I through III, 
the signal is shown by the thin black solid lines, the medians are shown by the thick black 
solid lines, and other quartiles are shown by the gray lines. Panel IV plots the differences 
between the third and the first quartiles of the outputs of the filters. The incoming signal 

10 (1/10 of the amplitude) is shown at the bottom of this panel. The distance between the 
time ticks is equal to the time constant of the filters T. 

29.2 Detection of Intermittency 

Fig. 56 provides an example of usage of AARFs for signal amplitudes and threshold cross- 
ings to detect intermittency. Panel I illustrates that outputs of AARFs for signal amplitudes 

15 and threshold crossing rates for a signal with intermittency can be substantially different. 
The quartile outputs (for q ~ 0.25, 0.5, and 0.75 quantiles) of an AARF for signal threshold 
crossing rates are shown by the solid black lines, and the respective outputs of an AARF 
for signal amplitudes, by dashed fines. Panel II shows the median outputs of AARFs for 
threshold crossing rates (black solid lines) and amplitudes (dashed lines), and Panel III 

20 plots the difference between these outputs. In Panels I and II, the input signal is shown by 
gray lines. 
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29.3 Removing Outliers (Filtering of Impulse Noise) 

One of the most appealing features of the rank filters is their insensitivity to outliers, 
although the definition of outliers is different for the accelerations, threshold crossings, and 
amplitudes. For signal amplitudes, the insensitivity to outliers means that sudden changes 

5 in the amplitudes of the signal x{t), regardless of the magnitude of these changes^ do not 
significantly affect the output of the filter Dq{t) unless these changes persist for about the 
(1 - g)th fraction of the width of the moving window. The example in Fig. 57 illustrates such 
insensitivity of median amplitude AARFs and ARSs to outliers. The original uncorrupted 
signal is shown by the thick black line in the upper panel, and the signal+noise total by 

10 a thinner line. In the middle panel, the noisy signal is filtered through an RCio Cauchy 
test function median AARF (thick line), and an averaging RCw filter with the same time 
constant (thinner line). The distance between the time ticks is equal to 10 T, where T is 
the time constant of the filters. In the lower panel, the signal is filtered through an ARS 
emulator of a 5-point digital median filter (thick line), and a 5-point running mean filter 

15 (thinner line). The distance between the time ticks is equal to 50 sampling intervals. 

Another example of insensitivity of the RCio moving window median filter of signal 
amplitudes to outlier noise is given in Fig. 58. Outlier noise (Panel I) is added to the signal 
shown in Panel II. The total power of the noise is more than 500 times larger than the 
power of the signal, but the noise affects only ^ 25% of the data points. The periodogram 

20 of the signal+noise total is shown in Panel III, and the periodogram of the signal only is 
shown in Panel IV. The composite signal is filtered through an ARS emulator of a 10-point 
digital median filter, and the periodogram of the result is shown in Panel V. 
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29.4 Comparison of Outputs of Digital Rank Order Filters with 
Respective Outputs of AARFs, ARSs, and ARFs Based on Ideal 
Measuring System 

As an additional example of the particular advantage of employing the "real" acquisition 
systems as opposed to the "ideal" systems, let us compare the output of a digital median 
filter with the respective outputs of AARF, ARS, and the median filter based on an ideal 
5 measuring system. 



29.4.1 Respective Numerical Algorithms 

A formal equation for an analog rank filter based on an ideal measuring system is obviously 
contained in the prior art. For example, the equation for a square window rank filter 
10 based on an ideal measuring system is derived as follows. The time-averaged output of the 
ideal discriminator is described by the function (Nikitin et al., 1998, p. 169, Eq. (38), for 
example) 

n{D, t) = {e[D~ x(5)])y = ^ ^* ^ds d[D- x{s)] , (146) 

which is formally a surface in a three-dimensional rectangular coordinate system, whose 
15 points with the coordinates D, and satisfy Eq. (146). Thus the output of a rank filter 
for the qi\x quantile is the level line S7(D,t) = g of this surface. Even though this surface 
is a discontinuous surface, it can be formally differentiated with respect to threshold using 
the relation d6{x)ldx = 6{x) (see, for example, Eq. (4)). Substitution of Eq. (146) into 
the equation for the level fine, Eq. (11), leads to the formal expression for the ideal rank 
20 filter in a square window as 

■ ^ dt njp, t) ^ e[D- x{t)] -eiD~x{t- t)] 

dDn{D,t) T {5[D-x{s)])^ ' ^ ^ 

which is practically useless since it does not contain q explicitly. This difficulty can be 
easily overcome by setting 
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^(^)]>t)a° 

which leads to the approximate expression for the partial time derivative of 0.{D, t) as 
1 



(148) 



dtn{D,t) ; 



AT 



{{9[D-x{s)])^~q} . 



(149) 



Since the denominator of Eq. (147) cannot be numerically computed, we replace 
{6 [D — x(s)])y by its unimodal approximation, for example, by the approximation of 
Eq. (70). The combination of this approximation with Eqs. (147) and (149) leads to the 
numerical algorithm for the square window (amplitude) rank filter, based on the ideal mea- 
suring system, which can read as follows: 



Dn-,1-- 



(150) 



fr.-h {(M - !)/.-! + ^ exp [-i^^] } 

where we set AT equal to MAt, At being the sampling interval. 

The respective numerical algorithms for the A ARE and ARS are 
based on the integration of Eqs. (117) and (129), respectively, by the 
forward Euler method (Press et al., 1992, Chapter 16, for example). 



15 

29.4.2 Pileup Signal 

As was discussed in Section 14 (see also the more detailed discussion in Nikitin, 1998, for 
example), the unimodal approximation of Eq. (70) should be adequate for a signal with 
strong pileup effects. Thus one would expect that a rank filter based on an ideal measuring 
20 system might be satisfactory for such a signal. Fig. 59 a compares the outputs of a digital 
median filter with the respective outputs of an AARF (Panel I), an ARS (Panel II), and an 
ARF based on an ideal measuring system (Panel III). In all panels, the pileup signal is shown 
by the gray lines, the outputs of the digital median filter are shown by the dashed black 
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lines, and the respective outputs of the analog median filters are shown by the sohd black 
lines. As can be seen from this figure, even though the AARF and the ARS outperform 
the median filter based on an ideal measuring system, the performance of the latter might 
still be considered satisfactory. 



5 

29.4.3 Asymmetric Square Wave Signal 

For the amplitude density of the asymmetric square wave signal, shown by the gray lines 
in Fig. 59 b, the unimodal approximation of Eq. (70) is a poor approximation. As a result, 
the median filter based on an ideal measuring system (solid black line in Panel III) fails to 

10 adequately follow the output of the digital median filter (dashed line in the same panel). 
Panels I and II of Fig. 59 b compare the outputs of the digital median filter with the 
respective outputs of an AARF (Panel I) and an ARS (Panel II). In all examples, AT was 
chosen as T/10, where T is the width of the rectangular window. This figure illustrates 
the advantage of employing the "real" acquisition systems over the solutions based on the 

15 "ideal" systems of the prior art. 



29.4.4 Comparison of RCiq AARF with RCiq ARF Based on Ideal Measuring System 

Since conventional digital rank order filters employ rectangular time windows, it is difficult 
to directly compare the outputs of such filters with the RCm adaptive analog rank filters. 
20 Panel I of Fig. 60 implements such comparison of the quartile outputs of a digital square 
window rank filter (dashed lines) with the respective outputs of the RCio AARF (solid black 
lines). Panel II of the same figure compares the quartile outputs of the digital rank filter 
(dashed lines) with the respective outputs of the RCw ARF, based on an ideal measuring 
system (soHd black lines). The time constants of the analog filters were chosen as T/2 
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in both examples. The incoming signal is shown by the gray lines. Since the amplitude 
density of this signal is bimodal, the unimodal approximation of Eq. (70) does not insure a 
meaningful approximation of the respective digital rank filter by the "ideal" ARF. 



30 Summary of Main Transformations 



30.1 Modulated Density and Cumulative Distribution 

Fig. la shows a simplified schematic of the basic AVATAR system, summarizing various 
transformations of an input variable into scalar field variables (e.g., into densities, cumula- 
tive distributions, or counting rates), such as the transformations described by Eqs. (52), 
(54), (55), (56), (60), (61), (62), (64), (131), and (138) in this disclosure. This system 
can be implemented through various physical means such as electrical or electro-optical 
hardware devices, as well as in computer codes (software). The detailed description of 
Fig. la is as follows: 

The system shown in Fig. la is operable to transform an input variable into an output 
variable having mathematical properties of a scalar field of the Displacement Variable. A 
Threshold Filter (a Discriminator or a Probe) is applied to a difference of the Displacement 
Variable and the input variable, producing the first scalar field of the Displacement Variable. 
This first scalar field is then filtered with a first Averaging Filter, producing the second 
scalar field of the Displacement Variable. Without optional modulation, this second scalar 
field is also the output variable of the system, and has a physical meaning of either an 
Amplitude Density (when the Threshold Filter is a Probe), or a Cumulative Amplitude 
Distribution (when the Threshold Filter is a Discriminator) of the input variable. 

A Modulating Variable can be used to modify the system in the following manner. 
First, the output of the Threshold Filter (that is, the first scalar field) can be multiplied 
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(modulated) by the Modulating Variable, and thus the first Averaging Filter is applied to 
the resulting modulated first scalar field. For example, when the Threshold Filter is a Probe 
and the Modulating Variable is a norm of the first time derivative of the input variable, 
the output variable has an interpretation of a Counting (or Threshold Crossing) Rate. The 

5 Modulating Variable can also be filtered with a second Averaging Filter having the same 
impulse response as the first Averaging Filter, and the output of the first Averaging Filter 
(that is, the second scalar field) can be divided (normalized) by the filtered Modulating 
Variable. The resulting output variable will then have a physical interpretation of either 
a Modulated Threshold Density (when the Threshold Filter is a Probe), or a Modulated 

10 Cumulative Threshold Distribution (when the Threshold Filter is a Discriminator). For 
example, when the Threshold Filter is a Probe and the Modulating Variable is a norm of 
the first time derivative of the input variable, the output variable will have an interpretation 
of a Counting (or Threshold Crossing) Density. 



30.2 Mean at Reference Threshold 

15 Fig. 61 illustrates such embodiment of AVATAR as the transformation of an input variable 
into a Mean at Reference Threshold variable (see Eqs. (53) and (140)). As has been 
previously described in this disclosure, a comparison of the Mean at Reference Threshold 
with the simple time (or spatial) average will indicate the interdependence of the input 
and the reference variables. This transformation can be implemented by various physical 

20 means such as electrical or electro- optical hardware devices, as well as in computer codes 
(software). The detailed description of Fig. 61 is as follows: 

The system shown in the figure is operable to transform an input variable into an 
output Mean at Reference Threshold variable. A Probe is applied to the difi'erence of 
the Displacement Variable and the reference variable, producing a first scalar field of the 

25 Displacement Variable. This first scalar field is then modulated by the input variable, 
producing a modulated first scalar field of the Displacement Variable. This modulated first 
scalar field is then filtered with a first Averaging Filter, producing a second scalar field of 
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the Displacement Variable. The first scalar field is also filtered with a second Averaging 
Filter having the same impulse response as the first Averaging Filter, and the output of the 
first Averaging Filter (that is, the second scalar field) is divided by the filtered first scalar 
field. The resulting quotient is the Mean at Reference Threshold variable. 



30.3 QuANTiLE Density, Quantile Domain Factor, and Quantile Volume 

Among various embodiments of AVATAR, the ability to measure (or compute from digital 
data) (1) Quantile Density, (2) Quantile Domain Factor, and (3) Quantile Volume for a 
variable are of particular importance for analysis of variables. Quantile Density indicates 
the value of the density likely to be exceeded, Quantile Domain contains the regions of the 
highest density, and Quantile Volume gives the (total) volume of the quantile domain. The 
definitions of these quantities and a means of their implementation are unavailable in the 
existing art. Fig. 62 provides a simplified schematic of transforming an input variable into 
output Quantile Density, Quantile Domain Factor, and Quantile Volume variables according 
to Eqs. (142), (144), and (145), Notice that these transformations can be implemented by 
various physical means such as electrical or electro-optical hardware devices, as well as in 
computer codes (software). The detailed description of Fig. 62 is as follows: 

The upper portion of Fig. 62 reproduces the system shown in Fig. la where the threshold 
filter is a Probe. The output of such system is either an Amplitude Density, or a Modulated 
Threshold Density (MTD). This density can be further transformed into Quantile Density, 
Quantile Domain Factor, and Quantile Volume as described below. 

A second Probe is applied to the difference between a feedback of the Quantile Den- 
sity variable and the Amplitude Density/MTD, producing a first function of the Quantile 
Density variable. This first function of Quantile Density is then multiplied by the Ampli- 
tude Density/MTD, producing a first modulated function of Quantile Density. The first 
modulated function of Quantile Density is then filtered with a first Time Averaging Filter 
producing a first time averaged modulated function of Quantile Density, and integrated over 
the values of the Displacement Variable producing a first threshold integrated function of 
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Quantile Density. 

A first Discriminator, which is a respective discriminator of the second Probe, is apphed 
to the difference between the feedback of the Quantile Density variable and the Amplitude 
Density/ MTD, producing a second function of the Quantile Density variable. A quantile 

5 value and the second function of Quantile Density is then subtracted from a unity, and 
the difference is multiplied by the AmpHtude Density/MTD. This produces a second mod- 
ulated function of Quantile Density. This second modulated function of Quantile Density 
is then filtered with a second Time Averaging Filter having the impulse response of the 
first derivative of the impulse response of the first Time Averaging Filter. This filtering 

10 produces a second time averaged modulated function of Quantile Density. This second 
time averaged modulated function is then integrated over the values of the Displacement 
Variable producing a second threshold integrated function of Quantile Density. By divid- 
ing the second threshold integrated function by the first threshold integrated function and 
time- integrating the quotient, the system outputs the Quantile Density variable. 

15 By applying a second Discriminator to the difference of the Amplitude Density/MTD 

and the Quantile Density variable, the latter variable is transformed into the Quantile 
Domain Factor variable. By integrating the Quantile Domain Factor over the values of the 
Displacement Variable, the system outputs the Quantile Volume variable. 



30.4 Rank Normalization 

20 Fig. 63 provides a simplified schematic of such important embodiment of AVATAR as 
Rank Normalization of an input variable with respect to a cumulative distribution function 
generated by a reference variable (see, for example, Eq. (86)). The system shown in Fig. 63 
can be implemented through various physical means such as electrical or electro-optical 
hardware devices, as well as in computer codes (software). The detailed description of 

25 Fig. 63 is as follows: 

A Discriminator is apphed to the difference of the Displacement Variable and the ref- 
erence variable producing a first scalar field of the Displacement Variable. This first scalar 
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field is then filtered with a first Averaging Filter, producing a second scalar field of the 
Displacement Variable. A Probe is applied to the difference of the Displacement Variable 
and the input variable producing a third scalar field of the Displacement Variable. This 
third scalar field is multiplied by the second scalar field and the product is integrated over 

5 the values of the Displacement Variable to output the Rank Normalized variable. 

A Modulating Variable can be used to modify the system as follows. First, the output of 
the Discriminator (that is, the first scalar field) is modulated by the Modulating Variable, 
and thus the first Averaging Filter is applied to the resulting modulated first scalar field. 
The Modulating Variable is also filtered with a second Averaging Filter having the same 

10 impulse response as the first Averaging Filter, and the output of the first Averaging Filter 
(that is, the second scalar field) is divided (normalized) by the filtered Modulating Variable. 
The resulting Rank Normalized variable will then have a physical interpretation of the input 
variable normalized with respect to a MTD of the reference variable. 



30.5 Explicit Analog Rank Filter 

15 Fig. 64 shows a schematic of an explicit Analog Rank Filter operable to transform an input 
scalar (or scalar field) variable into an output Rank Filtered variable according to Eq. (105) 
or Eq. (134). This filtering system can be implemented by various physical means such as 
electrical or electro-optical hardware devices, as well as in computer codes (software). The 
detailed description of Fig. 64 is as follows: 

20 A first Probe is applied to the difference of the Displacement Variable and the input 
variable producing a first scalar function of the Displacement Variable. This first scalar 
function is then filtered by a first Averaging Filter producing a first averaged scalar function 
of the Displacement Variable. A Discriminator, which is a respective discriminator of 
the first Probe, is applied to the difference of the Displacement Variable and the input 

25 variable, producing a second scalar function of the Displacement Variable. This second 
scalar function is then filtered with a second Averaging Filter having the same impulse 
response as the first Averaging Filter, producing a second averaged scalar function of the 
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Displacement Variable. 

A second Probe with a small width parameter is applied to the difference of a quantile 

value and the second averaged scalar function producing an output of the second Probe. 

This output is multiplied by the first averaged scalar function and by the Displacement 
5 Variable, This product is then integrated over the values of the Displacement Variable 

producing the output Rank Filtered variable. 

The first scalar function and the second scalar function can be also modulated by a 

Modulating Variable, and the first averaged scalar function and the second averaged scalar 

function can be divided by the Modulating Variable filtered with a third Averaging Filter, 
10 which third Averaging Filter has an impulse response identical to the impulse response of 

the first and second Averaging Filters. Then the Rank Filtered variable will correspond to a 

certain quantile of the Modulated Cumulative Threshold Distribution of the input variable. 



30.6 Feedback Analog Rank Filter 

Fig. 65 provides a simplified schematic of a feedback Analog Rank Filter for a single scalar 
15 variable or a scalar field variable, following the expressions of Eqs. (113) and (135). This 
filter can be embodied in various hardware devices such as electrical or electro-optical, or 
in computer codes (software). The detailed description of Fig. 65 is as follows: 

A Probe is applied to the difference between a feedback of the Rank Filtered variable and 
the input variable producing a first scalar function of the Rank Filtered variable. This first 
20 scalar function is filtered with a first Time Averaging Filter producing a first averaged scalar 
function of the Rank Filtered variable. A Discriminator, which is a respective discriminator 
of the Probe, is applied to the difference between the feedback of the Rank Filtered variable 
and the input variable producing a second scalar function of the Rank Filtered variable. 
This second scalar function is subtracted from a quantile value, and the difference is filtered 
25 with a second Time Averaging Filter having the impulse response of the first derivative of 
the impulse response of the first Time Averaging Filter, producing a second averaged scalar 
function of the Rank Filtered variable. The second averaged scalar function is divided by 
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the first averaged scalar function, and the quotient is time-integrated to output the Rank 
Filtered variable. 

The first scalar function, and the difference between the quantile value and the second 
scalar function, can also be modulated by a Modulating Variable. Then the Rank Filtered 

5 variable will correspond to a certain quantile of the Modulated Cumulative Threshold Dis- 
tribution of the input variable. 

The input variable can also be a scalar field variable. Then averaging by Spatial Aver- 
aging Filters with identical impulse responses may be added to the averaging by the first 
and second Time Averaging Filters. These Spatial Averaging Filters should be operable on 

10 the spatial coordinates of the input variable. When modulation by a Modulating Variable 
is implemented in a system for rank filtering of a scalar field variable, the Spatial Averaging 
Filters should be operable on the spatial coordinates of the input variable and on the spatial 
coordinates of the Modulating Variable. 



30.7 Analog Rank Filter for Ensemble of Scalar Variables 

15 Fig. 66 provides a diagram of a feedback Analog Rank Filter for a discrete ensemble or 
scalar variables, as described by Eq. (130). This filter can be embodied in various hardware 
devices such as electrical or electro-optical, or in computer codes (software). The detailed 
description of Fig. 66 is as follows: 

A Probe is applied to each difference between a feedback of the Rank Filtered variable 

20 and each component of the ensemble of input variables, producing a first ensemble of 
scalar functions of the Rank Filtered variable. Each component of the first ensemble of 
scalar functions is multiplied by the weight of the respective component of the ensemble 
of input variables, and the products are added together producing a first scalar function 
of the Rank Filtered variable. This first scalar function is then filtered with a first Time 

25 Averaging Filter, producing a first averaged scalar function of the Rank Filtered variable. 
A Discriminator, which is a respective discriminator of the Probe, is applied to each 
difference between the feedback of the Rank Filtered variable and each component of the 
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ensemble of input variables producing a second ensemble of scalar functions of the Rank 
Filtered variable. Each difference between a quantile value and each component of the 
second ensemble of scalar functions is multiplied by the weight of the respective component 
of the ensemble of input variables, and the products are summed, which produces a second 

5 scalar function of the Rank Filtered variable. This second scalar function is further filtered 
with a second Time Averaging Filter having the impulse response of the first derivative 
of the impulse response of the first Time Averaging Filter, producing a second averaged 
scalar function of the Rank Filtered variable. The second averaged scalar function is then 
divided by the first averaged scalar function, and the quotient is time-integrated to output 

10 the Rank Filtered variable. 

Optional modulation by an ensemble of Modulating Variables can be added to the 
system. Then the output of the Probe (that is, the first ensemble of scalar functions) 
is modulated by the ensemble of Modulating Variables (component-by-component), and 
the difference between the quantile value and the output of the Discriminator (that is, 

15 the difference between the quantile value and the second ensemble of scalar functions) is 
modulated by the ensemble of Modulating Variables (component- by-component). 
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BRIEF DESCRIPTION OF FIGURES 

The patent or application file contains at least one drawing executed in color. Copies of 
this patent or patent application publication with color drawing(s) will be provided by the 
Office upon request and payment of necessary fee. 

5 

Figure 1 a. Simplified schematic of basic system for analysis of variables. 
Figure 1 b. Simplified schematic of basic elements of system for analysis of variables. A 
scalar input variable x{t) (Panel I) is transformed by a discriminator (Panel Ila) and by a 
differential discriminator , or prohe (Panel lib), into continuous functions of two variables, 
10 displacement (threshold) D and time t, as shown in Panels Ilia and Illb. 

Figure 2. Input-output characteristics of some exemplary discriminators and the respec- 
tive probes (differential discriminators) . 

Figure 3. Illustration of the counting process for a continuous signal. The upper part 
of the figure shows a computer generated signal x{t) with crossings of the threshold D at 

15 times t^. The Heaviside step function of the difference of the signal x{t) and the threshold 
D is shown in the middle of the figure. The differential of the function d[x{t) - D] equals 
±1 at times and is shown at the bottom of the figure. Reproduced from (Nikitin, 1998). 
Figure 4. Introduction to Modulated Threshold Density. Consider intersections of 
a scalar variable (signal) x{t) in the interval [0,T] with the thresholds {Dj}, where 

20 Dj+i = Dj + AD. The instances of these crossings are labeled as {tj}, ti+i > ti. The 
thresholds {Dj} and the crossing times {ti} define a grid. We shall name a rectangle of 
this grid with the lower left coordinates (ti,Dj) as a Sij box. We will now identify the time 
interval At^j as tt+i — if the box s^^ covers the signal (as shown), and zero otherwise. 
Figure 5. Example of using modulated densities for measuring the input variable K in 

25 terms of the reference variable x. Notice that the amplitude densities of the fragments of 
the signals Xi{t) and x^it] shown in the left-hand panels are identical. Notice also that 
the modulating signals Ki{t), K2{2), and Ks{t) are identical for the respective modulated 
densities of the signals xi{t) and X2{t), while the modulated densities are clearly different. 
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Thus even though the amplitude densities and the modulating signals are identical, 
different reference signals still result in different modulated densities. 
Figure 6. Diagram illustrating an optical threshold smoothing filter (probe). 
Figure 7. Diagram illustrating transformation of an input variable into a modulated 
5 threshold density. 

Figure 8. RCin impulse response functions for n = 0 (exponential forgetting), n — l, and 
n = 2. 

Figure 9 a. Amplitude, counting, and acceleration densities of a signal. The left column 
of the panels shows the fragments of three different signals in rectangular windows. The 

10 second column of the panels shows the amplitude densities, the third column shows 
the counting densities, and the right column shows the acceleration densities for these 
fragments. This figure illustrates that the acceleration and counting densities generally 
reveal different features of the signal than do the amplitude densities. 
For the fragment xi{t) (the upper row of the panels), \x{t]\ = constant, and thus the 

15 counting and the amplitude densities are identical. For the fragment X2{t) (the middle row 
of the panels), \x{t) \ = constant, and thus the acceleration and the amplitude densities are 
identical. 

Figure 9 b. Amplitude, counting, and acceleration densities of a signal. The left column 
of the panels shows the fragments of three different signals in rectangular windows. The 

20 second column of the panels shows the amplitude densities, the third column shows 
the counting densities, and the right column shows the acceleration densities for these 
fragments. This figure illustrates that the acceleration and counting densities generally 
reveal different features of the signal than do the amplitude densities. 
Figure 10. Example of time dependent acceleration densities, threshold crossing rates, 

25 and amplitude densities computed in a 1-second rectangular moving window for two 
computer generated non-stationary signals (Panels la and lb). Panels 2a and 2b show the 
acceleration densities, Panels 3a and 3b show the threshold crossing rates, and Panels 4a 
and 4b show the amplitude densities. 

Figure 11. Illustration of applicability of quantile densities, domains, and volumes to 
30 analysis of scalar variables. 

Figure 12 a. Quantile densities, volumes, and domains displayed as time dependent 
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quantities computed in a 1-second rectangular sliding window for the signal shown in 
Panel la of Fig. 10. 

Figure 12 b. Quantile densities, volumes, and domains displayed as time dependent 
quantities computed in a 1-second rectangular sliding window for the signal shown in 

5 Panel lb of Fig. 10. 

Figure 13. Phase space densities of a signal. The first column of the patiels in the figure 
shows the fragments of three different signals in rectangular windows. The second column 
of the panels shows the phase space amplitude densities, and the third column displays the 
phase space counting densities. 

10 Figure 14. Example of time dependent phase space amplitude densities computed ac- 
cording to Eq. (60) in a 1-second rectangular moving window for two computer generated 
non-stationary signals shown in Panels la and lb of Fig. 10. The figure plots the level 
lines of the phase space amplitude densities (Panels la and 2a), at times indicated by the 
time ticks. Panels lb and 2b show the time slices of these densities at time t = to- 

15 Figure 15. Example of time dependent phase space counting rates computed according 
to Eq. (62) in a 1-second rectangular moving window for two computer generated non- 
stationary signals shown in Panels la and lb of Fig. 10. The figure plots the level lines of 
the phase space counting rates (Panels la and 2a), at times indicated by the time ticks. 
Panels lb and 2b show the time slices of these rates at time t = to- 

20 Figure 16. Boundaries of the median domains for the phase space amplitude densities. 
The upper panel shows the boundary for the signal of Panel la of Fig. 10, and the lower 
panel shows the median domain boundary for the signal of Panel lb of Fig. 10. 
Figure 17. Boundaries of the median domains for the phase space counting densities. 
The upper panel shows the boundary for the signal of Panel la of Fig. 10, and the lower 

25 panel shows the median domain boundary for the signal of Panel lb of Fig. 10. 
Figure 18. Schematic statement of the underlying motivation behind AVATAR. 
Figure 19. Simplified conceptual schematic of a possible hardware device for displaying 
time dependent amplitude densities of a scalar variable. 

Figure 20. Simplified conceptual schematic of a possible hardware device for displaying 
30 time dependent threshold crossing rates of a scalar variable. 

Figure 21. Illustration of possible hardware device for displaying time slices of phase 
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space amplitude densities. 

Figure 22. Illustration of possible hardware device for displaying time slices of phase 
space counting rates. 

Figure 23. Estimator H^(D, t) of Eq. (63) in ^ = 9/10 quantile domain, computed for the 
5 two computer generated nonstationary scalar signals shown in Panels la and lb. Panels 2a 
and 2b display the values of the estimator for K = \x\, and Panels 3a and 3b display these 
values for K = 

Figure 24. Illustration of adequacy of the approximation of Eq, (73) when the signals 
x{t) and y{t) represent responses of linear detector systems to trains of pulses with high 

10 incoming rates, Poisson distributed in time. 

Figure 25, Illustration of the resulting density as a convolution of the component densities 
for uncorrelated signals. The signals Xi{t), X2{t), and Xi{t) + X2{t) are shown in the left 
column of the panels, and the respective panels in the right column show the respective 
amplitude densities. The signal X2{t} is random (non-Gaussian) noise. In the lower right 

15 panel, the measured density of the combined signal is shown by the solid line, and the 
density computed as the convolution of the densities bi{D) and 62 (-D) is shown by the 
dashed line. 

Figure 26. The amplitude {5 {D - x))^ and the counting {{\x\)j)'^ {\x\5 {D - x))j, 
densities of the fragment of the signal shown in the upper panel. One can see that the 
20 Gaussian unimodal approximation (dashed lines) is more suitable for the counting density 
than for the amplitude density. 

Figure 27. Example of the usage of the estimator given by Eq, (92) for quantification of 
changes in a signal. The signals are shown in Panels la and lb. The distributions Ca{D, t) 
are computed in a 1-second rectangular moving window as the amplitude (for Panels 2a 

25 and 2b) and counting (for Panels 3a and 3b) cumulative distributions. Thus yg{t) are 
the outputs of the respective rank filters for these distributions. The estimators Qabi^', q] 
are computed as the outputs of the Gaussian normalizer of Eq, (83), The values of these 
outputs for different quartile values are plotted by the gray (for q = 1/2), black (for 
q = 1/4), and light gray (for q — 3/4), In this example, the estimator Qab(t;q) quantifies 

30 the deviations of Ca(D, t) from the respective normal distributions. 

Figure 28. Example of usage of rank normalization for discriminating between different 



117 



Atty. Docket No.' 15262.2 
First Named Inventor: Nikitm 
US Express: EL388608922US 

pulse shapes of a variable. Panel I shows the input signal consisting of three different 
stretches, 1 through 3, corresponding to the variables shown in Fig. 9. Panel 11 displays 
the difference between C|^°,,(x,t) and Ci°{x,t), where the reference signal r is a Gaussian 
process with the mean Kw and the variance K20 — KIq, and K^m a^re computed for the 

5 input signal x{t). This difference is zero for the first stretch of the input signal, since for 
this stretch the amplitude and the counting densities are identical (see Fig. 9). Panel III 
displays the difference between C^^i^^{x,t) and Ci°{x,t). This difference is zero for the 
second stretch of the input signal, since for this stretch the amplitude and the acceleration 
densities are identical {see Fig. 9). The distance between the time ticks is equal to the 

10 constant T of the time filter. 

Figure 29. Additional example of sensitivity of the difference between two rank normal- 
ized signals to the nature of the reference distributions. Panel I shows the input signal, 
and Panel II displays the amplitude density computed in an RCio window with the time 
constant T, equal to the distance between the time ticks. Panel III plots the difference 

15 between C^^^^^{x,t) and Ci°{x,t), where the reference signal r is a Gaussian process with 
the mean Kio and the variance if 20 — K^w, computed for the input signal x{t). 
Figure 30. Simplified flowchart of an analog rank normalizer. 

Figure 31. Simplified flowchart of a device for comparison of two signals by a means of 
rank normalization. 

20 Figure 32. Time averages of the absolute values of the differences, {\C'l^.^{x,t) - 
Ci°x{x,t)\)T, for K = \x\ and K = \x\, for the variable shown above the panel of the 
figure. The distance between the time ticks is equal to the time constants T of the filtering 
windows. 

Figure 33. Illustration of sensitivity of the amplitude and counting phase space densities 
25 to differences in the signal's wave form. The panels in the left column show the sound 
signals for several letters of the alphabet. The top signals in the individual panels are the 
original input signals. The normahzed input signals and their normalized first derivatives, 
respectively, are plotted below the original input signals. The middle column of the 
panels shows the amplitude, and the right column the counting densities of these pairs of 
30 normalized signals. 

Figure 34. Panel I: The original speech signal "Phase Space" is shown in the top of 
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the panel. This signal is normalized with respect to a Gaussian process with the mean 
and variance of the original signal in a moving rectangular window of 45 ms, and the 
result is plotted just below the original signal. The bottom of the panel shows the time 
derivative of the speech signal, normalized the same way. Panel 11: Time slices of the 

5 threshold density c{Dx-, Dy,t), where x and y are the normalized original signal and its 
normalized derivative, respectively, and c{D^, Dy,t) is their amplitude density in the time 
window 45 ms. The slices are taken approximately through the middles of the phonemes. 
Panel III: Time slices of the cumulative distribution C{Dx, Dy,t), where x and y are the 
normahzed original signal and its normalized derivative, respectively, and C{Dx,Dy,t) is 

10 their distribution in the time window 45 ms. The slices are taken approximately through 
the middles of the phonemes. Panel IV: The value of the estimator of a type of Eq. (94), 
where the reference distribution is taken as the average distribution computed in the 
neighborhood of the phonemes "a" . 

Figure 35. Outline of an optical speech recognition device. 

15 Figure 36. Illustration of the relationship between the outputs of a rank filter and the level 
fines of the amplitude distribution of a scalar signal. Panel I shows the input signal x{t) 
on the time-threshold plane. This signal can be viewed as represented by its instantaneous 
density 5[D —x[t)]. Threshold integration by the discriminator Tad{D) transforms this 
instantaneous density into the threshold averaged distribution T/^d[D — ^if)] (Panel II). 

20 This distribution is further averaged with respect to time, and the resulting distribution 
B{D,t) — {J^AD [D — x{s)])j, is shown in Panel III. The quartile level lines are computed 
as the outputs of the rank filter given by Eq. (105), and are plotted in the same panel. 
Panel IV shows the input signal x{t), the level lines of the amplitude distribution for 
q = 1/4, 1/2, and 3/4 (gray lines), and the outputs of a digital rank order filter (black 

25 lines). 

Figure 37. Example of Fig. 36, repeated for the respective analog and digital median 
filters for the discrete input signals. The instantaneous density of a discrete signal can be 
represented by 5[D - x{t)] J2iS{t - U), as shown in Panel I. Panel II shows the threshold 
averaged distribution J^ad[D - x{t)]Y.^5{t - ti), and Panel III of the figure compares 
30 the level line B{D,t) = {Tad [D - x{s)] Y.i5{s - ti))j, = 1/2 (solid black line) with the 
respective output of a digital median filter (white dots). 
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Figure 38. Simplified schematic of a device for analog rank filtering. 
Figure 39. Simplified schematic of a device for analog rank filtering. Module I of the 
device outputs the signal (iT)^ [g - C'}^{Dq,t)], and Module II estimates {K)l^c),{D^,t). 
The outputs of Modules I and II are divided to form Dq{t), which is integrated to produce 
the output of the filter Dq{t). 

Figure 40. Simplified schematic of the implementation of Eq. (113) in an analog device, 
with the adaptation according to Eq. (111). Module I takes the outputs of Modules II 
and III as inputs. The output of Module I is also a feedback input of Module II. Module IV 
outputs ^D{t), which is used as one of the inputs of Module II (the width parameter of 
the discriminator and the probe) for adaptation. 

Figure 41. Comparison of the quartile outputs (for q = 0.25, 0.5, and 0.75 quantiles) 
of the Cauchy test function RCn AARF for signal amplitudes with the corresponding 
conventional square window digital order statistic filter. The outputs of the AARF are 
shown by the thick black solid lines, and the respective outputs of the square window 
order statistic filter are shown by the thin black lines. The time constant of the impulse 
response of the analog filter is T, and the corresponding width of the rectangular window 
is 2aT, where a is the solution of the equation a - ln(l + a) = ln(2). The incoming signal 
is shown by the gray line, and the distance between the time ticks is equal to 2aT. 
Figure 42. Comparison of the quartile outputs (for q = 0.25, 0.5, and 0.75 quantiles) of 
the Cauchy test function square window AARF for signal amplitudes with the correspond- 
ing conventional square window digital order statistic filter. The outputs of the AARF 
are shown by the black solid lines, and the respective outputs of the square window order 
statistic filter are shown by the dashed lines. The widths of the time windows are T in all 
cases. The incoming signal is shown by the gray line, and the distance between the time 
ticks is equal to T. 

Figure 43. Finding a rank of a discrete set of numbers according to Eq. (124). Five 
numbers Xi are indicated by the dots on the X-axis of the top panel. The solid line shows 
the density resulting from the spatial averaging with a Gaussian test function, and the 
dashed lines indicate the contributions into this density by the individual numbers. The 
solid fine in the middle panel plots the cumulative distribution. The crosses indicate rCg(a) 
and :F^D[xg(a)] at the successive integer values of the parameter a. The bottom panel 
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plots the evolution of the value of Xq{a) in relation to the values of Xi. 
Figure 44. Analog rank selection for an ensemble of variables. In Panel I, the solid line 
shows the 3rd octile of a set of four variables (xi{t) through 0:4 (t), dashed Hnes), computed 
according to Eq. (129). In Panel II, the solid line shows the median {q = 1/2 in Eq. (129)) 
5 of the ensemble. The thick dashed line plots the median digitally computed at each 
sampling time. The time constant of the analog rank selector is ten times the sampling 
interval. 

Figure 45. Comparison of the quartile outputs (for q = 0.25, 0.5, and 0.75 quantiles) of a 
square window digital order statistic filter (dashed lines) with its emulation by the Cauchy 
10 test function ARS (solid black lines). The incoming signal is shown by the gray line, and 
the distance between the time ticks is equal to the width of the time window T. 
Figure 46. Simplified schematic of a device (according to Eq. (129)) for analog rank 
selector for three input variables. 

Figure 47. Example of performance of AARFs for ensembles of variables. This figure also 
15 illustrates the fact that counting densities do not only reveal different features of the signal 
than do the amplitude densities, but also respond to different changes in the signal. The 
figure shows the outputs of median AARFs for an ensemble of three variables. The input 
variables are shown by the gray lines. The thicker black lines in Panels I and II show the 
outputs of the median AARFs for amplitudes, and the thinner black lines in both panels 
20 show the outputs of the median AARFs for counting densities. All AARFs employ Cauchy 
test function and RCw time averaging. The distance between the time ticks in both panels 
is equal to the time constant of the time filters. 

Figure 48. Diagram illustrating transformation of a scalar field into a modulated 
threshold density. 

25 Figure 49. Simplified schematic of a device (according to Eq. (133)) for analog rank filter 
of a discrete monochrome surface with 3x3 spatial averaging. 

Figure 50. Filtering out static impulse noise from an image according to the algorithm of 
Eq. (137). Panel 1: The original image Z. Panel 2: The image corrupted by a random 
unipolar impulse noise of high magnitude. About 50% of the image is affected. Panel 3a: 
30 The initial condition for the filtered image is a plane of constant magnitude. Panels 3b 
through 3g: The snapshots of the filtered image Q (the first decile of the corrupted one, 
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q = 1/10) at steps n. 

Figure 51. Filtering out time- varying impulse noise according to the algorithm of 
Eq. (137). Panels la through Ic: Three consecutive frames of an image corrupted 
by a random (bipolar) impulse noise of high magnitude. About 40% of the image 
is affected. Panels 2a through 2c: The image filtered through a smoothing filter, 
= 'Lm,nW.mnZi^rn,3~n- Panels 3a through 3c: The rank filtered image Q (the 
median, q = 1/2). The smoothing filter in Eq. (137) is the same used in Panels 2a 
through 2c. 

Figure 52. Diagram illustrating transformation of a vector field into a modulated 
threshold density. 

Figure 53. Diagram of a process for the transformation of the incoming vector field 
x(a,i) into a modulated threshold density c;^(D;a,i), and the subsequent evaluation of 
the quantile density Zq{t), quantile domain factor 5'^(D;a,i), and the quantile volume 
-^^(a, t) of this density. 

Figure 54 a. Comparison of the median densities and volumes computed directly from 
the definitions (Eqs. (24) and (25), gray lines) with those computed through Eqs. (143) 
and (145) (black lines). Panels la and 2a relate to the amplitude densities, and Panels 3a 
and 4a relate to the counting densities. 

Figure 54 b. Comparison of the median densities and volumes computed directly from 
the definitions (Eqs. (24) and (25), gray lines) with those computed through Eqs. (143) 
and (145) (black fines). Panels lb and 2b relate to the amplitude densities, and Panels 3b 
and 4b relate to the counting densities. 

Figure 55. Quartile outputs (for q = 1/4 through 3/4 quantiles) of the RCi^ Cauchy test 
function AARFs for the signal amplitudes (Panel I), threshold crossing rates (Panel II), and 
threshold crossing accelerations (Panel III). The signal consists of three different stretches, 
1 through 3, corresponding to the signals shown in Fig. 9. In Panels I through III, the 
signal is shown by the thin black solid lines, the medians are shown by the thick black 
solid lines, and other quartiles are shown by the gray lines. Panel IV plots the differences 
between the third and the first quartiles of the outputs of the filters. The incoming signal 
(1/10 of the amplitude) is shown at the bottom of this panel. The distance between the 
time ticks is equal to the time constant of the filters T. 
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Figure 56. Detection of intermittency. Panel I illustrates that outputs of AARFs for 
signal amplitudes and threshold crossing rates for a signal with intermittency can be 
substantially different. The quartile outputs (for q = 0.25, 0.5, and 0.75 quantiles) of 
an AARF for signal threshold crossing rates are shown by the solid black hnes, and the 
5 respective outputs of an AARF for signal amplitudes, by dashed lines. Panel 11 shows the 
median outputs of AARFs for threshold crossing rates (black solid lines) and amplitudes 
(dashed lines), and Panel III plots the difference between these outputs. In Panels I and II, 
the input signal is shown by gray lines. 

Figure 57. Insensitivity of median amplitude AARFs and ARSs to outliers. The original 
10 uncorrupted signal is shown by the thick black line in the upper panel, and the signal+noise 
total by a thinner line. In the middle panel, the noisy signal is filtered through an RCio 
Cauchy test function median AARF (thick line), and an averaging RCio filter with the 
same time constant (thinner line). The distance between the time ticks is equal to 10 T, 
where T is the time constant of the filters. In the lower panel, the signal is filtered through 
15 an ARS emulator of a 5-point digital median filter (thick line), and a 5-point running mean 
filter (thinner line). The distance between the time ticks is equal to 50 sampling intervals. 
Figure 58. Outlier noise (Panel I) is added to the signal shown in Panel II. The total 
power of the noise is more than 500 times larger than the power of the signal, but the 
noise affects only 25% of the data points. The periodogram of the signal+noise total 
20 is shown in Panel III, and the periodogram of the signal only is shown in Panel IV. The 
composite signal is filtered through an ARS emulator of a 10-point digital median filter, 
and the periodogram of the result is shown in Panel V. 

Figure 59 a. Comparison of the outputs of a digital median filter (dashed lines) with 
the respective outputs of an AARF, an ARS, and an ARFs based on an ideal measuring 

25 system (solid lines), for a signal (gray lines) with strong pileup effects. 

Figure 59 b. Comparison of the outputs of a digital median filter (dashed lines) with 
the respective outputs of an AARF, an ARS, and an ARFs based on an ideal measuring 
system (solid lines), for an asymmetric square wave signal (gray lines). 
Figure 60. Comparison of the quartile outputs of a digital square window rank filter 

30 (dashed lines in both panels) with the respective outputs of the RCio AARF (solid black 
fines in Panel I), and with the quartile outputs of the RCiq ARF, based on an ideal 
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measuring system (solid black lines in Panel II) . 

Figure 61. Schematic of transforming an input variable into an output Mean at Reference 
Threshold variable. 

Figure 62. Schematic of transforming an input variable into output Quantile Density, 
5 Quantile Domain Factor, and Quantile Volume variables. 

Figure 63. Schematic of Rank Normalization of an input variable with respect to a 
reference variable. 

Figure 64. Schematic of an explicit Analog Rank Filter. 

Figure 65. Schematic of an Analog Rank Filter for a single scalar variable or a scalar 
10 field variable. 

Figure 66. Schematic of an Analog Rank Filter for an ensemble of scalar variables. 
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BEST MODES FOR CARRYING OUT THE INVENTION 
1 MAIN EQUATIONS FOR PRACTICAL EMBODIMENTS 

(i) Modulated threshold density, Eq. (52): 

(ii) Mean at reference threshold, Eq. (53): 

( /r [D - X(S)] )y 

(iii) Amplitude density, Eq. (54): 

6(D, t)=^{fl ^A^AD, [A - Xi{s)] 

(iv) Counting density, Eq. (55): 



(v) Counting rates, Eq. (56): 



(vi) Phase space amplitude density, Eq. (60): 
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(vii) Phase space counting density, Eq. (61): 

( + ^o^^^o^ - dD.T^D^ [Di: - x{s)] ) 
r{D,,Di,t) = ^ LT 



(viii) Phase space counting rates, Eq. (62): 

7^(£>„ Di, t) = {^^ {x Dif + [x D^f Qd^T^v^ [D^ - x{s)] dv.J'^D, [Di - x{s)]^ . 

(ix) Estimator of differences in quantile domain between the mean at reference thresh- 
old and the time average, Eq. (63): 

S,(D.)^ l™g-;>-<'^)^l .[</.(D-x)).-/,Wl . 

(x) Modulated cumulative distribution, Eq. (64): 

^-oo {K{s))j, 

(xi) Rank normalization with respect to the reference distribution Ck,v{^: t), Eq. (86): 

(xii) Rank normalization by a discriminator with an arbitrary input-output response, 
Eq. (88): 

y{t) = ^R.(t) [Dr{t) - X(£)] = n^AD,„W [Dr,i{t) - Xi{t)] . 

(xiii) Rank normalization of a scalar variable by a discriminator with an arbitrary 
input-output response, Eq. (89): 

=^ ^^ 2iK,o-Kl) l^ioit) - x{t)] . 

(xiv) Estimator of differences between two distributions, Eq. (92): 

Qat{t;q) = C,[ygit),t] 



(xv) Statistic for comparison of two distributions, Eq. (95): 

Kbit) = /_Vx/l(x)i?[Ca(x,t),C'&(x,i),Ca(x,t),C5(x,t)] . 

(xvi) Statistic of Cramer- von Mises type, Eq. (97): 

Kbit) = y°°dC„(x,^)u;[Ca(x,^)] W[Ca{x,t) - Cb{x,t)] = 

= y"^d''xca(x,i) w[Ca{x,t)] W[Cai^,t) - Cb{x,t)] . 

(xvii) Probabilistic comparison of amplitudes, Eq. (102): 

m = Ijisg^it - s) Cl% = (|cff. J . 

(xviii) Analog rank filter, Eq. (105): 

D,it) = l2dDDcKiD,t)S[CKiD,t)-q]^ 
« {DcKiD,t)d,J'AACKiD,t)-q])Z . 

(xix) Adaptive analog rank filter, Eq. (113): 

^ _ q {Kis)f^ - {Kjs) J^^njs) - xis)])l 

{Kis)dnJ'AVis)ms)-xis)]Y^ 

(xx) Alternative embodiment of adaptive analog rank filter, Eq. (117): 

. _ qiKis))!^-{Kis):F^nis ) [D,it) - xis)])^ 

U — ■ — r— . 



AT 



[{Kis)dnJ'c.Dis) [D,it)-xis)])\ 



(xxi) Threshold averaged instantaneous density for a continuous ensemble of variables, 
Eq. (120): 

6(D;i,n(/^)) =y^%n(/i)/R[D-x^(t)] . 
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(xxii) Threshold averaged instantaneous cumulative distribution for a continuous en- 
semble of variables, Eq. (121): 

BiD;t,n(fi)) = |_^d/.n(/i) [D - x^(i)] . 

(xxiii) Modulated density for a continuous ensemble of variables, Eq. (122): 

C;,{D;t,nM) = jj^n(^) j^^-^j^ E . 

(xxiv) Modulated cumulative distribution for a continuous ensemble of variables, 
Eq. (123): 

C^(D;t,nM)=y_^d^„M X . 

(xxv) Analog rank selector for a continuous ensemble, Eq. (126): 

• Jf^'^dg dtQTjt, e) { BnM [a, )^ 

(xxvi) RCio analog rank selector for a discrete ensemble, Eq. (129): 

q-Ezni J^AD [Xqjt) - Xi(t)] 



T(Ein,ao^Ai>K(s)-x,(s)])J° ' 

(xxvii) Adaptive analog rank filter for a discrete ensemble of variables, Eq. (130): 

(E. K,{s) dnJ^ABis) ms) - x,{s)])l 

(xxviii) Modulated threshold density for a scalar field, Eq. (131): 

, „ _ {K{':.^)dnrAD [P~^(r, s)] )»i 
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(xxix) RCio analog rank selector /filter for a scalar field (n-dimensional surface), 
Eq. (133): 

i _ g-(-^AD[^,(x,t)-^(r,t)])4 

T{{dnJ'^nM:^,s)~z{r,s)])i)l° ' 

(xxx) Analog rank filter for a scalar field, Eq. (134): 

D,(a, t) = l^dD D ck{D- a, t) d^J^^g [^^(D; a, t) - q] . 

(xxxi) Adaptive analog rank filter for a scalar field, Eq. (135): 

. , , 9 {K{v, s)) - {K{r, s) J^AD(.,s) s) - xir, s)]f/^ 

ID 13. — ~— f I ^J\. 

(xxxii) Numerical algorithm for analog rank processing of an image given by the matrix 
Z = Z,j(t), Eq. (136): 

Qk = Qfc-i + {q- F)/ffc 

F = T:m,nn^n.nJ'AD [Qk-l - {Z,-m,j -n) k-^l] 
g = Em,n «^m« ^iP^AD [Qfc-1 - (^i-m ,i-n)fe-l] 

(xxxiii) Modulated threshold density for a vector field, Eq. (138): 

c.(D;M)=^^^-"^^-^^--^-)J>^"-. 

(i^(r,.))J;i 

(xxxiv) Modulated threshold density for an ensemble of vector fields, Eq. (139): 

^(ir,(r,.)/R[D-x,(r,g)])^;^^ 



CK(D;3L,t,n{iJ,)) = J djin{ij)- 



(xxxv) Mean at reference threshold for a vector field input variable, Eq. (140): 
(K(r,s)/R[D-x(r,5)])^,, 



{MxK}T,A(D;a,i) = 



(/R[D-x(r,s)])^_^ 
( n^i Ki(r, s) dn^J^AD. [Dj - Xjjr, s)] 
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(xxxvi) Analog quantile density filter, Eq. (142): 

{f^^dr^Tz{r,s)dn:F^D[zg{s)-z{r,s)])''^ 

(xxxvii) Analog quantile domain filter, Eq. (144): 

5,(D; a, t) - J^^n [z{B, t) - Zg(t)] . 

(xxxviii) Analog quantile volume filter, Eq. (145): 

R,{ai,t)= l2^d''TSgir;a,t) = (5,(r;a,i))^ . 



2 Articles of Manufacture 

Various embodiments of AVATAR may include hardware, firmware, and software embod- 
iments, that is, may be wholly constructed with hardware components, programmed into 
firmware, or be implemented in the form of a computer program code. 

Still further, the invention disclosed herein may take the form of an article of manu- 
facture. For example, such an article of manufacture can be a computer-usable medium 
containing a computer-readable code which causes a computer to execute the inventive 
method. 



130 



